Online adaptive discrete empirical interpolation for nonlinear model - - PowerPoint PPT Presentation

Online adaptive discrete empirical interpolation for nonlinear model reduction

Benjamin Peherstorfer (Courant Institute, NYU) and Karen Willcox (MIT)

1 / 22

Introduction: Nonlinear PDEs with parameters

Discretized equations stemming from (steady-state) nonlinear PDE A

• N×N

y(µ) + f (y(µ))

• f :RN→RN

= 0

◮ Spatial domain Ω ⊂ R{1,2,3} ◮ Linear term A ∈ RN×N ◮ State variable y ∈ RN ◮ Parameter µ ∈ D ⊂ Rd ◮ Nonlinear term f (y(µ)) ∈ RN

2 / 22

Introduction: Discrete Empirical Interpolation Method

Derive POD-Galerkin reduced system

◮ Collect snapshots

y(µ1), . . . , y(µM) ∈ RN

◮ Construct POD basis V ∈ RN×n, n ≪ N of reduced space ◮ Project onto reduced space

V TAV

n×n

˜ y(µ) + V T

• n×N

f ( V

• N×n

˜ y(µ)) = 0 DEIM interpolates f as linear combination of basis U

◮ Compute nonlinear snapshots f (y(µ1)), . . . , f (y(µM)) ∈ RN ◮ Compute DEIM basis U ∈ RN×m of nonlinear snapshots ◮ Select interpolation points in P ∈ RN×m at which to sample f

V TAV

n×n

˜ y(µ) + V TU(PTU)−1

• n×m

PTf (V ˜ y(µ)) = 0

[Barrault et al., 2004], [Grepl et al., 2007], [Astrid et al., 2008], [Chaturantabut et al., 2010], [Carlberg et al., 2011], [Drohmann et al., 2012], [Drmac, Gugercin, 2016] 3 / 22

Nonlinear approximation

Have t ∈ [1, 3], grid x = [x1, . . . , x64]T in [0, 2π] f (t) = t (sin(xt) + sin(xtπ/2) + sin(xtπ)) Dimension of manifold induced by f is 1 M = {f (t) | t ∈ [1, 3]} ⊂ R64 Classical reduced models approximate M with space U = span{u1, . . . , um} Here, one-dimensional space U = span{u1} fails to approximate M ⇒ Localized and adaptive reduced spaces

4 / 22

2 4 6 −5 5 spatial domain solution time t = 1 function f 2 4 6 −5 5 spatial domain mode DEIM basis vector 2 4 6 −5 5 spatial domain approximation DEIM approx

Nonlinear approximation

Have t ∈ [1, 3], grid x = [x1, . . . , x64]T in [0, 2π] f (t) = t (sin(xt) + sin(xtπ/2) + sin(xtπ)) Dimension of manifold induced by f is 1 M = {f (t) | t ∈ [1, 3]} ⊂ R64 Classical reduced models approximate M with space U = span{u1, . . . , um} Here, one-dimensional space U = span{u1} fails to approximate M ⇒ Localized and adaptive reduced spaces

4 / 22

2 4 6 −5 5 spatial domain solution time t = 1 function f 2 4 6 −5 5 spatial domain mode DEIM basis vector 2 4 6 −5 5 spatial domain approximation DEIM approx

Nonlinear approximation

Have t ∈ [1, 3], grid x = [x1, . . . , x64]T in [0, 2π] f (t) = t (sin(xt) + sin(xtπ/2) + sin(xtπ)) Dimension of manifold induced by f is 1 M = {f (t) | t ∈ [1, 3]} ⊂ R64 Classical reduced models approximate M with space U = span{u1, . . . , um} Here, one-dimensional space U = span{u1} fails to approximate M ⇒ Localized and adaptive reduced spaces

4 / 22

2 4 6 −5 5 spatial domain solution time t = 2.98 2 4 6 −5 5 spatial domain mode 2 4 6 −5 5 spatial domain approximation function f DEIM basis vector DEIM approx

from global reduced spaces to localized subspaces U ⇒ U1, . . . , Uk

◮ Eftang, Patera, Ronquist; 2010 ◮ Dihlmann, Drohmann, Haasdonk; 2011 ◮ Amsallem, Zahr, Farhat; 2012 ◮ Eftang, Stamm; 2012 ◮ Eftang, Patera; 2013 ◮ Maday, Stamm; 2013 ◮ P., Butnaru, Willcox, Bungartz; 2014 ◮ Zimmermann; 2014 ◮ Albrecht, Haasdonk, Kaulmann, Ohlberger; 2015 ◮ Schindler, Ohlberger; 2015, 2017 ◮ Buhr, Engwer, Ohlberger, Rave; 2017a,2017b ◮ ...

5 / 22

reduced models with online adaptive (DEIM) subspaces U ⇒ U1, . . . , Uk ⇒ U(t)

◮ Amsallem, Zahr, Washabaugh; 2015 ◮ Carlberg; 2015 ◮ Kramer, P., Willcox; 2017 ◮ Lass; 2014 ◮ P., Willcox; 2015a, 2015b, 2016 ◮ Schindler, Ohlberger; 2015, 2017 ◮ Zahr, Farhat; 2015 ◮ ...

6 / 22

reduced models with online adaptive (DEIM) subspaces U ⇒ U1, . . . , Uk ⇒ U(t)

◮ Amsallem, Zahr, Washabaugh; 2015 ◮ Carlberg; 2015 ◮ Kramer, P., Willcox; 2017 ◮ Lass; 2014 ◮ P., Willcox; 2015a, 2015b, 2016 ◮ Schindler, Ohlberger; 2015, 2017 ◮ Zahr, Farhat; 2015 ◮ ... P., Willcox: Online Adaptive Model Reduction for Nonlinear Systems via Low-Rank Updates. SIAM Journal on Scientific Computing, 37(4):A2123-A2150, SIAM, 2015. ◮ Sparse residual from full model ⇒ avoid pre-computed quantities ◮ Low-rank updates ⇒ online efficiency ◮ Additive updates ⇒ arbitrary changes to reduced basis

6 / 22

Nonlinear approximation through online adaptivity

... ... . reduced state vector ˜ y(µM+i) sample f samples of non- linear function process DEIM interpolant approximate nonlinear function approximation

• uter loop

7 / 22

Nonlinear approximation through online adaptivity

• uter loop iteration

... ... . reduced state vector ˜ y(µM+1) sample f samples of non- linear function process DEIM interpolant approximate nonlinear function approximation

7 / 22

Nonlinear approximation through online adaptivity

• uter loop iteration

... ... . reduced state vector ˜ y(µM+1) sample f samples of non- linear function process DEIM interpolant approximate nonlinear function approximation reduced state vector ˜ y(µM+2) sample f samples of non- linear function process approximate nonlinear function approximation

7 / 22

Nonlinear approximation through online adaptivity

• uter loop iteration

... ... . reduced state vector ˜ y(µM+1) sample f samples of non- linear function process DEIM interpolant approximate nonlinear function approximation reduced state vector ˜ y(µM+2) sample f samples of non- linear function process adapted DEIM interpolant approximate nonlinear function approximation adapt

7 / 22

Nonlinear approximation through online adaptivity

• uter loop iteration

... ... . reduced state vector ˜ y(µM+1) sample f samples of non- linear function process DEIM interpolant approximate nonlinear function approximation reduced state vector ˜ y(µM+2) sample f samples of non- linear function process adapted DEIM interpolant approximate nonlinear function approximation adapt reduced state vector ˜ y(µM+3) sample f samples of non- linear function process adapted DEIM interpolant approximate nonlinear function approximation adapt ...

7 / 22

Nonlinear approximation through online adaptivity

• uter loop iteration

... ... . reduced state vector ˜ y(µM+1) sample f samples of non- linear function process DEIM interpolant approximate nonlinear function approximation reduced state vector ˜ y(µM+2) sample f samples of non- linear function process adapted DEIM interpolant approximate nonlinear function approximation adapt reduced state vector ˜ y(µM+3) sample f samples of non- linear function process adapted DEIM interpolant approximate nonlinear function approximation adapt ...

Key ingredients of online adaptation

◮ Sparse residual from full model ⇒ avoid pre-computed quantities ◮ Low-rank updates ⇒ online efficiency ◮ Additive updates ⇒ arbitrary changes to DEIM space

7 / 22

Adaptive DEIM: Problem setup

Outer loop

◮ Online phase consists of k = 1, . . . , M′ outer loop iterations ◮ Reduced solution ˜

y(µM+k) requested at each iteration

◮ DEIM interpolant used to approximate f (V ˜

y(µM+k)) Adaptive interpolant

◮ DEIM interpolant (U0, P0) built offline from snapshots ◮ Adaptation initialized with (U0, P0) ◮ In each step k = 1, . . . , M′

◮ Adapt DEIM basis Uk−1 to obtain Uk ◮ Adapt DEIM interpolation points matrix Pk−1 to obtain Pk ◮ Derive adapted interpolant (Uk, Pk) ◮ Use (Uk, Pk) to approximate nonlinear term 8 / 22

Adaptive DEIM: Problem setup

Outer loop

◮ Online phase consists of k = 1, . . . , M′ outer loop iterations ◮ Reduced solution ˜

y(µM+k) requested at each iteration

◮ DEIM interpolant used to approximate f (V ˜

y(µM+k)) Adaptive interpolant

◮ DEIM interpolant (U0, P0) built offline from snapshots ◮ Adaptation initialized with (U0, P0) ◮ In each step k = 1, . . . , M′

◮ Adapt DEIM basis Uk−1 to obtain Uk ◮ Adapt DEIM interpolation points matrix Pk−1 to obtain Pk ◮ Derive adapted interpolant (Uk, Pk) ◮ Use (Uk, Pk) to approximate nonlinear term 8 / 22

Adaptive DEIM: Oversampling

DEIM interpolates f at interpolation points {p1, . . . , pm}, i.e.,

• PT

U(PTU)−1PTf (y(µ)) − f (y(µ))

• 2 = 0

Oversample nonlinear function

◮ Draw ms ∈ N points uniformly from {1, . . . , N} \ {p1, . . . , pm} ◮ First m points equal DEIM interpolation points ◮ Create sampling points matrix S ∈ RN×m+ms ◮ Solve regression problem to obtain coefficient

c(y(µ)) = (STU)+STf (y(µ))

◮ Residual

r(y(µ)) = Uc(y(µ)) − f (y(µ)) In general, non-zero residual at sampling points, STr(y(µ))2 > 0

9 / 22

Adaptive DEIM: Low-rank basis updates

Window contains previous w ∈ N reduced solutions ˜ y(µk1), . . . , ˜ y(µkw ) Adapt space Uk−1 ∈ RN×m with rank-one update αkβT

k ∈ RN×m

Uk = Uk−1 + αkβT

k

Update αkβT

k minimizes

ST

k ((Uk−1 + αkβT k )C k − F k)2 F ◮ Sk ∈ RN×m+ms is sampling points matrix at step k ◮ C k = [c(V ˜

y(µk1)), . . . , c(V ˜ y(µkw ))] ∈ Rm×w coefficient matrix

◮ F k = [f (V ˜

y(µk1)), . . . , f (V ˜ y(µkw ))] ∈ RN×w RHS matrix With Rk = Uk−1C k − F k obtain minimization problem arg min

αk∈RN,βk∈Rm ST k Rk m+ms×w

+ ST

k αk m+ms

βT

k C k w

2

F

10 / 22

Adaptive DEIM: Construction of basis update

Optimization problem With Rk = Uk−1C k − F k obtain arg min

αk∈RN,βk∈Rm ST k Rk m+ms×w

+ ST

k αk m+ms

βT

k C k w

2

F

Update αkβT

k derived from m × m generalized eigenproblem ◮ Size m × m scales with dimension m of DEIM space ◮ Symmetric positive definite ◮ Optimal update that minimizes residual in · F

Costly steps of optimization

◮ Sampling nonlinear function at m + ms × w components ◮ QR transformation (necessary if C k ∈ Rm×w rank deficient) ◮ Generalized eigenproblem of size m × m

11 / 22

Adaptive DEIM: Interpolation points update

Point adaptation

◮ Offline point selection too expensive online ◮ Often not necessary to recompute all points ◮ Heuristic selects point to be replaced

Online interpolation points update

◮ Find basis vector that was rotated most

diag(UT

k Uk−1) ◮ Replace the corresponding interpolation point ◮ Rerun offline point selection for this point only ◮ Costs in O(Nm

• select

+ m3

• new

)

12 / 22

Adaptive DEIM: Toy example

13 / 22

2 4 6 −5 5 spatial domain solution time t = 1 target solution 2 4 6 −5 5 spatial domain mode static basis vector 2 4 6 −5 5 spatial domain approximation static DEIM 2 4 6 −5 5 spatial domain adapted basis vector 2 4 6 −5 5 spatial domain adaptive DEIM

Adaptive DEIM: Toy example

13 / 22

2 4 6 −5 5 spatial domain solution time t = 1 target solution 2 4 6 −5 5 spatial domain mode static basis vector 2 4 6 −5 5 spatial domain approximation static DEIM 2 4 6 −5 5 spatial domain adapted basis vector 2 4 6 −5 5 spatial domain adaptive DEIM

Numerical results: Optimization problem

Define Ω = D = [0, 1]2 and g : Ω × D → R g(x, µ) = µ1µ2 exp(x1x2) exp(20x − µ2

2) ◮ Find parameter µ ∈ D that maximizes

• Ω

g(x, µ)dx

◮ Discretize in Ω on 40 × 40 equidistant grid ◮ Approximate nonlinear function with DEIM ◮ Search for optimum with DEIM interpolant 0.2 0.4 0.6 0.8 1 µ1 0.2 0.4 0.6 0.8 1 µ2 50 100 150 200 250 value of objective

1e-06 1e-04 1e-02 1e+00 4 8 12 16 20 L2 error of maxima #DEIM basis vectors built offline

14 / 22

Numerical results: Optimization problem

1e-06 1e-05 1e-04 1e-03 1e-02 1e-01 1e+00 100 200 300 400 500

• ptimization error (14)

adaptivity step k adapt, m = 5 basis vectors static, m = 100 static, m = 5 1e-05 1e-04 1e-03 120 130 140 150 160 170

• ptimization error (14)
• nline time [s]

static, 65-100 DEIM basis vec. adapt, 200-300 samples (a) error (b) speedup ◮ Search for maximum with Nelder-Mead algorithm ◮ DEIM interpolant adapted when Nelder-Mead evaluates function ◮ Adaptive space with 5 dim achieves accuracy of 100-dim static space ◮ Adaptive interpolant reduces online runtime

15 / 22

Numerical results: FitzHugh-Nagumo system

Models electrical activity in a neuron ǫ∂ty v = ǫ2∂2

x y v + f (y v) − y w + c ,

∂ty w = by v − γy w + c ,

0.2 0.4 0.6 0.8 1 time 0.25 0.5 0.75 1 x

• 0.4
• 0.2

0.2 nonlinear function

◮ Spatial domain x ∈ Ω = [0, 1], time domain t ∈ T = [0, 1] ◮ Voltage y v : Ω × T → R and recovery of voltage y w : Ω × T → R ◮ Nonlinear function f (y v) = y v(y v − 0.1)(1 − y v) ◮ Equidistant grid in Ω, N = 2048 DoFs ◮ Forward Euler in time domain, 106 time steps

[Chaturantabut et al., 2010] 16 / 22

Numerical results: FitzHugh-Nagumo system

1e-06 1e-04 1e-02 2 4 6 8 10 12 rel L2 error, averaged over time #DEIM basis vectors static adapt, 100 samples adapt, 200 samples 1e-06 1e-04 1e-02 56 58 60 62 64 rel L2 error, averaged over time

• nline time [s]

static adapt, 100 samples adapt, 200 samples (a) error (b) speedup ◮ Adapt DEIM interpolant every 200-th time step ◮ Accuracy improvement of up to two orders of magnitude ◮ Overall runtime not dominated by evaluation of nonlinear term ◮ Therefore minor speedup compared to static DEIM interpolant

17 / 22

Numerical results: Combustor model

Nonlinear advection-diffusion-reaction κ∆y − ν∇y + f (y, µ) = 0 in Ω

◮ f (µ, t) = [y H2, y O2, y H2O, T]T with

◮ mass fractions y H2, y O2, y H2O ◮ temperature T

◮ Nonlinear term f with two parameters

◮ pre-exponential factor µ1 = A ◮ activation energy µ2 = E x1 [mm] x2 [mm] 0.5 1 1.5 0.2 0.4 0.6 0.8 temp [K] 500 1000 1500 2000

µ = (5.5 · 1011, 1.5 · 103)

x1 [mm] x2 [mm] 0.5 1 1.5 0.2 0.4 0.6 0.8 temp [K] 500 1000 1500 2000

µ = (1.5 · 1013, 9.5 · 103)

[Buffoni & Willcox, 2010] 18 / 22

Numerical results: Region of interest

1e-06 1e-05 1e-04 1e-03 1000 2000 3000 4000 5000 rel L2 error in region of interest adaptivity step k ms = 40 samples ms = 60 samples ms = 80 samples 1e-06 1e-05 1e-04 1e-03 1e-02 15000 16000 17000 18000 19000 20000 rel L2 error in region of interest

• nline time [s]

static adapt (a) error (b) speedup ◮ Adapt interpolant to region of interest DRoI ⊂ D ◮ Accuracy improvement of up to three orders of magnitude ◮ Significant online runtime improvement through adaptation

19 / 22

Numerical results: Extended/shifted parameter domain

2e+03 1e+04 5.5e+11 1.5e+13 parameter µ2 parameter µ1

• ffline domain

2e+03 1e+04 5.5e+11 1.5e+13 parameter µ2 parameter µ1

• ffline domain

(a) extended (b) shifted ◮ Adapted interpolant to extended/shifted parameter domain ◮ Static DEIM interpolant fails to approximate nonlinear term ◮ Accuracy improvement of up to three orders of magnitude ◮ Similar accuracy as rebuilding interpolant from scratch

20 / 22

Numerical results: Extended/shifted parameter domain

1e-07 1e-06 1e-05 1e-04 D0

E

D2

E

D4

E

D6

E

D8

E

D10

E

rel L2 error in extended domain domain static rebuilt adapt

2e+03 1e+04 5.5e+11 1.5e+13 parameter µ2 parameter µ1

• ffline domain

(a) extended (b) shifted ◮ Adapted interpolant to extended/shifted parameter domain ◮ Static DEIM interpolant fails to approximate nonlinear term ◮ Accuracy improvement of up to three orders of magnitude ◮ Similar accuracy as rebuilding interpolant from scratch

20 / 22

Numerical results: Extended/shifted parameter domain

1e-07 1e-06 1e-05 1e-04 D0

E

D2

E

D4

E

D6

E

D8

E

D10

E

rel L2 error in extended domain domain static rebuilt adapt 1e-07 1e-06 1e-05 1e-04 1e-03 D0

S

D2

S

D4

S

D6

S

D8

S

D10

S

rel L2 error in shifted domain domain static rebuilt adapt (a) extended (b) shifted ◮ Adapted interpolant to extended/shifted parameter domain ◮ Static DEIM interpolant fails to approximate nonlinear term ◮ Accuracy improvement of up to three orders of magnitude ◮ Similar accuracy as rebuilding interpolant from scratch

20 / 22

Numerical results: Expected failure rate

1e-04 1e-03 1e-02 1e+02 1e+03 1e+04 1e+05 RMSE #samples for mean static adapt static, IS adapt, IS 1e-04 1e-03 1e-02 1e+02 1e+03 1e+04 1e+05 RMSE #samples for mean static adapt static, IS adapt, IS (a) four DEIM basis vectors (b) six DEIM basis vectors ◮ Computed expected failure rate of combustor (temp above 2240K) ◮ Adapted interpolant while Monte Carlo sampling proceeds ◮ RMSE about an order of magnitude lower

21 / 22

Conclusions

Adaptive discrete empirical interpolation method

◮ Adaptivity targets nonlinear term approximation ◮ Avoids pre-computed quantities that restrict update (cf. localization) ◮ Adaptation is computationally cheap (sparse samples, low-rank updates)

Online adaptive model reduction

◮ Nonlinear approximation of solution manifold ◮ More robust with respect to initial reduced model construction ◮ Identifies problem structure that is amenable to low-rank approximation

References

◮

P., Willcox: Online Adaptive Model Reduction for Nonlinear Systems via Low-Rank Updates. SIAM Journal on Scientific Computing, 37(4):A2123-A2150, SIAM, 2015.

◮

P., Butnaru, Willcox, Bungartz: Localized Discrete Empirical Interpolation Method. SIAM Journal on Scientific Computing, 36(1):A168-A192, SIAM, 2014. 22 / 22