The Loewner Framework for Model Reduction of Flow Equations - - PowerPoint PPT Presentation

the loewner framework for model reduction of flow
SMART_READER_LITE
LIVE PREVIEW

The Loewner Framework for Model Reduction of Flow Equations - - PowerPoint PPT Presentation

The Loewner Framework for Model Reduction of Flow Equations Matthias Heinkenschloss Department of Computational and Applied Mathematics Rice University, Houston, Texas heinken@rice.edu February 18, 2020 Workshop Mathematics of Reduced Order


slide-1
SLIDE 1

The Loewner Framework for Model Reduction of Flow Equations

Matthias Heinkenschloss

Department of Computational and Applied Mathematics Rice University, Houston, Texas heinken@rice.edu

February 18, 2020 Workshop Mathematics of Reduced Order Models ICERM, Brown University, February 17 - 21, 2020 Joint work with A. C. Antoulas (ECE, Rice U. & MPI Magdeburg) and I. V. Gosea (MPI Magdeburg) Funded in part by NSF grants CCF-1816219 and DMS-1819144

Matthias Heinkenschloss

  • Feb. 18, 2020

1

slide-2
SLIDE 2

Motivation

◮ Full order model (FOM) (typically discretization of system of PDEs) E d dty(t) = Ay(t) + F(y(t)) + Bu(t) with state y(t) ∈ Rn, n ≫ 1, input u(t) ∈ Rm, and output z(t) = Cy(t) + Du(t) ∈ Rp. ◮ Projection based reduced order model (ROM): Matrices V, W ∈ Rn×r, r ≪ n; ROM WT EV d dt y(t) = WT AV y(t) + WT F(V y(t)) + WT Bu(t) with state y(t) ∈ Rr, r ≪ n, input u(t) ∈ Rm, and output

  • z(t) = CV

y(t) + Du(t) ∈ Rp. ◮ Goal: Input-to-output map u → z of ROM approximates input-to-output map u → z of FOM. ◮ Notation: States y, input u, output z.

Matthias Heinkenschloss

  • Feb. 18, 2020

2

slide-3
SLIDE 3

Many ROM methods, incl. ◮ Proper Orthogonal Decomposition (POD) (e.g., books and survey articles [Gubisch and Volkwein, 2017], [Hesthaven et al., 2015], [Quarteroni et al., 2016]). ◮ Reduced Basis (RB) methods (e.g., books and survey articles [Haasdonk, 2017], [Hesthaven et al., 2015], [Quarteroni et al., 2016], [Rozza et al., 2008]). ◮ Balanced Truncation and Balanced POD (e.g., book [Antoulas, 2005] and survey articles [Benner and Breiten, 2017], [Rowley and Dawson, 2017]). ◮ Rational Interpolation (e.g., book [Antoulas et al., 2020a] and survey article [Beattie and Gugercin, 2017]). All require access to E, A, . . . to generate projection matrices V, W ∈ Rn×r, then generate ROM E = WT EV, A = WT AV, . . . This talk: Loewner framework ◮ Constructs ROM E, A, . . . directly from data. ◮ For LTI systems [Antoulas et al., 2020a], [Antoulas et al., 2017]; recent work on descriptor systems, incl. Oseen equations (e.g., [Antoulas et al., 2020b], [Gosea et al., 2020]). ◮ Extension to quadratic-bilinear systems [Gosea and Antoulas, 2018], Burger’s equation [Antoulas et al., 2018].

Matthias Heinkenschloss

  • Feb. 18, 2020

3

slide-4
SLIDE 4

Loewner for Linear Time-Invariant (LTI) Systems

◮ Consider linear time-invariant system E d dty(t) = Ay(t) + b u(t) with state y(t) ∈ Rn, n ≫ 1, input u(t) ∈ R, and output z(t) = cT y(t) + d u(t) ∈ R. Example: Oseen equations. ◮ To simplify presentation consider case of single input b ∈ Rn, u(t) ∈ R, and single output c ∈ Rn, z(t) ∈ R (SISO). Multiple inputs, multiple outputs (MIMO) can be handled via so-called left and right tangential directions, but is a bit more technical and requires more notation.

Matthias Heinkenschloss

  • Feb. 18, 2020

4

slide-5
SLIDE 5

Frequency Domain

◮ Frequency domain s ∈ iR s E y(s) = Ay(s) + b u(s), z(s) = cT y(s). From now on we work in frequency domain; use same notation y, u, z for variables in frequency domain. ◮ Transfer function H(s) = cT (sE − A)−1b. ◮ Seek ROM E = WT EV, A = WT AV, b = WT b, c = VT c, s E y(s) = A y(s) + b u(s),

  • z(s) =

cT y(s) + d u(s), so that transfer function

  • H(s) =

cT (s E − A)−1 b ≈ H(s) = cT (sE − A)−1b.

Matthias Heinkenschloss

  • Feb. 18, 2020

5

slide-6
SLIDE 6

Interpolation

◮ Given distinct frequencies µ1, . . . µr, λ1, . . . λr ∈ C, want ROM s.t.

  • cT (µj

E − A)−1 b = H(µj) = H(µj) = cT (µjE − A)−1b, j = 1, . . . r,

  • cT (λj

E − A)−1 b = H(λj) = H(λj) = cT (λjE − A)−1b, j = 1, . . . r. ◮ Theorem: Interpolation guaranteed if ROM projection matrices V, W satisfy (λjE − A)−1b ∈ R(V), j = 1, . . . r,

  • cT (µjE − A)−1T ∈ R(W),

j = 1, . . . r. ◮ Can choose (assume these have full rank) V =

  • (λ1E − A)−1b, . . . , (λrE − A)−1b
  • ∈ Cn×r,

W∗ =    cT (µ1E − A)−1 . . . cT (µrE − A)−1    ∈ Cr×n.

Matthias Heinkenschloss

  • Feb. 18, 2020

6

slide-7
SLIDE 7

Elementary identities

Recall H(s) = cT (sE − A)−1b. cT (µE − A)−1E (λE − A)−1b = 1 µ − λ

  • cT (µE − A)−1(µ − λ)E (λE − A)−1b
  • =

1 µ − λ

  • cT (µE − A)−1

(µE − A) − (λE − A)

  • (λE − A)−1b
  • =

1 µ − λ

  • H(λ) − H(µ)
  • .

cT (µE − A)−1A (λE − A)−1b = −cT (µE − A)−1(µE − A) (λE − A)−1b + µcT (µE − A)−1E (λE − A)−1b = −H(λ) + µcT (µE − A)−1E (λE − A)−1b = −H(λ) + µ 1 µ − λ

  • H(λ) − H(µ)
  • =

1 µ − λ

  • λH(λ) − µH(µ)
  • .

Matthias Heinkenschloss

  • Feb. 18, 2020

7

slide-8
SLIDE 8

If V =

  • (λ1E − A)−1b, . . . , (λrE − A)−1b
  • ∈ Cn×r,

W∗ =    cT (µ1E − A)−1 . . . cT (µrE − A)−1    ∈ Cr×n, then (recall H(s) = cT (sE − A)−1b)

  • E = W∗EV = −

   

H(µ1)−H(λ1) µ1−λ1

· · ·

H(µ1)−H(λr) µ1−λr

. . . ... . . .

H(µr)−H(λ1) µr−λ1

· · ·

H(µr)−H(λr) µr−λr

   

def

= − L,

  • A = W∗AV = −

   

µ1H(µ1)−H(λ1)λ1 µ1−λ1

· · ·

µ1H(µ1)−H(λr)λr µ1−λr

. . . ... . . .

µrH(µr)−H(λ1)λ1 µr−λ1

· · ·

µrH(µr)−H(λr)λr µr−λr

   

def

= − Ls,

  • b = W∗b = [H(µ1), . . . , H(µr)]T ,
  • c = V∗c = [H(λ1), . . . , H(λr)]T .

Matthias Heinkenschloss

  • Feb. 18, 2020

8

slide-9
SLIDE 9

If V =

  • (λ1E − A)−1b, . . . , (λrE − A)−1b
  • ∈ Cn×r,

W∗ =    cT (µ1E − A)−1 . . . cT (µrE − A)−1    ∈ Cr×n, then (recall H(s) = cT (sE − A)−1b)

  • E = W∗EV = −

   

H(µ1)−H(λ1) µ1−λ1

· · ·

H(µ1)−H(λr) µ1−λr

. . . ... . . .

H(µr)−H(λ1) µr−λ1

· · ·

H(µr)−H(λr) µr−λr

   

def

= − L,

  • A = W∗AV = −

   

µ1H(µ1)−H(λ1)λ1 µ1−λ1

· · ·

µ1H(µ1)−H(λr)λr µ1−λr

. . . ... . . .

µrH(µr)−H(λ1)λ1 µr−λ1

· · ·

µrH(µr)−H(λr)λr µr−λr

   

def

= − Ls,

  • b = W∗b = [H(µ1), . . . , H(µr)]T ,
  • c = V∗c = [H(λ1), . . . , H(λr)]T .

Only need data (µj, H(µj)), (λj, H(λj)), j = 1, . . . , r, to construct ROM E, A, b, c!

Matthias Heinkenschloss

  • Feb. 18, 2020

8

slide-10
SLIDE 10

Loewner Framework

◮ Given distinct frequencies µ1, . . . µk, λ1, . . . λk ∈ C and transfer function measurements H(µ1), . . . , H(µk), H(λ1), . . . , H(λk) ∈ C, ◮ use Loewner matrix L =    

H(µ1)−H(λ1) µ1−λ1

· · ·

H(µ1)−H(λk) µ1−λk

. . . ... . . .

H(µk)−H(λ1) µk−λ1

· · ·

H(µk)−H(λk) µk−λk

    ∈ Ck×k and shifted Loewner matrix Ls =    

µ1H(µ1)−H(λ1)λ1 µ1−λ1

· · ·

µ1H(µ1)−H(λk)λk µ1−λk

. . . ... . . .

µkH(µk)−H(λ1)λ1 µk−λ1

· · ·

µkH(µk)−H(λk)λk µk−λk

    ∈ Ck×k ◮ to directly compute ROM E, A ∈ Rr×r, b, c ∈ Rr, such that

  • cT (µj

E − A)−1 b = H(µj) ≈ H(µj) = cT (µjE − A)−1b, j = 1, . . . k,

  • cT (λj

E − A)−1 b = H(λj) ≈ H(λj) = cT (λjE − A)−1b, j = 1, . . . k.

Matthias Heinkenschloss

  • Feb. 18, 2020

9

slide-11
SLIDE 11

Loewner Framework

◮ Given distinct frequencies µ1, . . . µk, λ1, . . . λk ∈ C and transfer function measurements H(µ1), . . . , H(µk), H(λ1), . . . , H(λk) ∈ C, ◮ use Loewner matrix L =    

H(µ1)−H(λ1) µ1−λ1

· · ·

H(µ1)−H(λk) µ1−λk

. . . ... . . .

H(µk)−H(λ1) µk−λ1

· · ·

H(µk)−H(λk) µk−λk

    ∈ Ck×k and shifted Loewner matrix Ls =    

µ1H(µ1)−H(λ1)λ1 µ1−λ1

· · ·

µ1H(µ1)−H(λk)λk µ1−λk

. . . ... . . .

µkH(µk)−H(λ1)λ1 µk−λ1

· · ·

µkH(µk)−H(λk)λk µk−λk

    ∈ Ck×k ◮ to directly compute ROM E, A ∈ Rr×r, b, c ∈ Rr, such that

  • cT (µj

E − A)−1 b = H(µj) ≈ H(µj) = cT (µjE − A)−1b, j = 1, . . . k,

  • cT (λj

E − A)−1 b = H(λj) ≈ H(λj) = cT (λjE − A)−1b, j = 1, . . . k. ◮ Typically more data k than ROM size r.

Matthias Heinkenschloss

  • Feb. 18, 2020

9

slide-12
SLIDE 12

Loewner Framework: Redundant Data

◮ Given distinct frequencies µ1, . . . µk, λ1, . . . λk ∈ C and transfer function measurements H(µ1), . . . , H(µk), H(λ1), . . . , H(λk) ∈ C. ◮ Compute Loewner L ∈ Ck×k and shifted Loewner Ls ∈ Ck×k. ◮ Compute (short) SVDs [L Ls] = Y1Σ1X∗

1,

  • L

Ls

  • = Y2Σ2X∗

2,

where Σ1 ∈ Rk×2k, Σ2 ∈ R2k×k, Y1, X2 ∈ Ck×k. ◮ Y, X ∈ Ck×r obtained by selecting first r columns of Y1 and X2. (r determined based on singular values.) ◮ ROM:

  • E = −Y∗LX,
  • A = −Y∗LsX,
  • b = Y∗

H(µ1), . . . , H(µk) T ,

  • cT =
  • H(λ1), . . . , H(λk)
  • X.

Matthias Heinkenschloss

  • Feb. 18, 2020

10

slide-13
SLIDE 13

Loewner Framework: Complex → Real

◮ Data µj, λj, . . . complex ⇒ previous ROM E, A, b, c complex. ◮ If data also includes conjugate complex data, {µj}k

j=1 = {µj}k j=1,

{λj}k

j=1 = {λj}k j=1,

then can transform complex ROM into real ROM with same transfer function.

Matthias Heinkenschloss

  • Feb. 18, 2020

11

slide-14
SLIDE 14

Loewner Framework: Transfer Function at Infinity

◮ So far considered FOMs s E y(s) = Ay(s) + b u(s), z(s) = cT y(s) with transfer function H(s) = Hspr(s) = cT (sE − A)−1b ◮ However, in general transfer function composed of strictly proper part and polynomial part H(s) = Hspr(s) + Hpoly(s), where Hspr(s) → 0 as |s| → ∞, and Hpoly(s) = d1 + d1s + . . .. In particular discretized Oseen has Hpoly(s) = 0. ◮ In practice, H(s) ≈ H(s) in range of frequency data, not outside. Loewner ROM generates transfer function

  • H(s) =

Hspr(s) ( Hpoly(s) = 0). ◮ Can compute Hpoly(s) or estimate (Loewner with high freq. data), and apply Lowener with measurements Hspr(s) = H(s) − Hpoly(s), to generate a ROM with transfer function

  • H(s) =

H (s) + H (s).

Matthias Heinkenschloss

  • Feb. 18, 2020

12

slide-15
SLIDE 15

Oseen Equations

∂ ∂tv(x, t) + (a(x) · ∇)v(x, t) − ν∆v(x, t) + ∇p(x, t) = 0 in Ω × (0, T), ∇ · v(x, t) = 0 in Ω × (0, T), (−p(x, t)I + ν∇v(x, t))n(x) = 0

  • n Γn × (0, T),

v(x, t) = 0

  • n Γd × (0, T),

v(x, t) = g(x, t)

  • n Γg × (0, T),

v(x, 0) = 0 in Ω, with ν = 1/50 and domain

1 2 3 4 5 6 7 8 0.5 1

and 6 inputs (suction blowing at inflow, backward facing step). Output:(weak form approximation of) z(t) =

  • Γobs
  • − p(x, t)I + ν∇v(x, t)
  • n(x)ds

( 2 outputs).

Matthias Heinkenschloss

  • Feb. 18, 2020

13

slide-16
SLIDE 16

Taylor-Hood FEM Semidiscretization of Oseen Equations

E11 d dtv(t) = A11v(t) + A12p(t) + B10g(t) + B11 d dtg(t), t ∈ (0, T), 0 = AT

12v(t) + B20g(t),

t ∈ (0, T), v(0) = 0, z(t) = C1v(t) + C2p(t) + D0g(t) + D1 d dtg(t), t ∈ (0, T).

d dtg(t) terms arise, e.g., when Dirichlet control inputs are applied.

Transfer function H(s) = Hspr(s) + P0 + s P1, with

P0 =D0 − C1E−1

11 A12(AT 12E−1 11 A12)−1B2,0 − C2(AT 12E−1 11 A12)−1AT 12E−1 11 B3,

P1 =D1 − C2(AT

12E−1 11 A12)−1

AT

12E−1 11 B1,1 + B2,0

  • ,

where B3 := B1,0 − A11E−1

11 A12(AT 12E−1 11 A12)−1B2,0.

Matthias Heinkenschloss

  • Feb. 18, 2020

14

slide-17
SLIDE 17

◮ Transfer function H(s) = Hspr(s) + P0 + s P1, has polynomials part induced by algebraic constraint

P0 =D0 − C1E−1

11 A12(AT 12E−1 11 A12)−1B2,0 − C2(AT 12E−1 11 A12)−1AT 12E−1 11 B3,

P1 =D1 − C2(AT

12E−1 11 A12)−1

AT

12E−1 11 B1,1 + B2,0

  • ,

where B3 := B1,0 − A11E−1

11 A12(AT 12E−1 11 A12)−1B2,0.

◮ Tedious to compute. ◮ Can estimate P0 and P1 from high frequency data. ◮ Instead of applying Loewner to H(s) = Hspr(s) + P0 + s P1 apply Loewner to (estimated) strictly proper part H(s) − P0 + s P1.

Matthias Heinkenschloss

  • Feb. 18, 2020

15

slide-18
SLIDE 18

Estimation of P0 and P1 in H(s) = Hspr(s) + P0 + s P1

◮ From one frequency η ≫ 1 H(ı η) = Hspr(ı η)

  • →0 (η→∞)

+P0 + ı ηP1 ≈ P0 + ı ηP1. Estimates P0 = ℜ

  • H(ı η)
  • and

P1 = η−1ℑ

  • H(ı η)
  • .

◮ From two frequencies θ > η ≫ 1. H(ı θ) − H(ı η) =

  • Hspr(ı θ) + P0 + ı θP1
  • Hspr(ı η) + P0 + ı ηP1
  • = Hspr(ı θ) − Hspr(ı η) + (ı θ − ı η)P1 ≈ (ı θ − ı η)P1.

Estimate

  • P1 = ℜ

H(ı θ) − H(ı η) ı θ − ı η

  • ... and

P0 = ℜ ı θH(ı θ) − ı ηH(ı η) ı θ − ı η

  • .

◮ Arbitrary number of (high) frequencies. Form Loewner Lhi and shifted Loewner Lhi

s . Estimate

  • P0 = ℜ
  • Lhi†Lhi

s

  • Rhi†

,

  • P1 = ℜ
  • Lhi†Lhi

Rhi† , (Lhi , Rhi matrices of directions.)

Matthias Heinkenschloss

  • Feb. 18, 2020

16

slide-19
SLIDE 19

Full Order Model # velocities nv = 12, 504 # pressures np = 1, 669 # inputs m = 6 # outputs p = 2 Loewner ROM Data: 20 logarithmically spaced points be- tween 103 and 105 and their complex conju- gates (will be able to construct real ROMs) to estimate P0 and P1. 100 logarithmically spaced points between 10−2 and 101 and their complex conjugates to compute ROM

20 40 60 80 100 120 140 160 180 200 10-15 10-10 10-5 100

Normalized singular val- ues of Loewner matrices [L Ls] and L Ls

  • .

ROM size r = 24 is small- est r with σr/σ1 > 10−10.

Matthias Heinkenschloss

  • Feb. 18, 2020

17

slide-20
SLIDE 20

Input 1 - Output 1

10-2 100 102 104 106 10-1 100 101 102 103 Original model Loewner - r = 24 Loewner new - r = 24 10-2 100 102 104 106 10-10 10-5 100

Loewner Loewner new

Input 3 - Output 2

10-2 100 102 104 106 100 105 Original model Loewner - r = 24 Loewner new - r = 24 10-2 100 102 104 106 10-10 10-5 100

Loewner Loewner new

Input 5 - Output 2

10-2 100 102 104 106 10-2 100 102 104 Original model Loewner - r = 24 Loewner new - r = 24 10-2 100 102 104 106 10-10 10-5 100

Loewner Loewner new

Matthias Heinkenschloss

  • Feb. 18, 2020

18

slide-21
SLIDE 21

Loewner for Quadratic-Bilinear Systems

Consider quadratic-bilinear SISO system E d

dty(t) = Ay(t) + G(y(t), y(t)) + Ny(t)u(t) + bu(t),

t ∈ (0, T), z(t) = cT y(t) + du(t), t ∈ (0, T), y(0) = 0, E, A, N ∈ Rn×n, b, c ∈ Rn, d ∈ R, and G : Rn × Rn → Rn is bilinear. System called quadratic-bilinear because y → G(y, y) is quadratic and (y, u) → Nyu is bilinear. Again present SISO case, but can extend to MIMO. Idea for ROM construction ◮ Expand quadratic-bilinear system into series of linear systems (e.g., [Rugh, 1981]). ◮ Identify transfer functions from expanded linear systems. ◮ Extend ROM methods for LTI system to generate ROM that approximates these transfer functions. ◮ Other approaches, e.g., [Benner and Breiten, 2015], [Ahmad et al., 2017]. ◮ Loewner [Antoulas et al., 2018].

Matthias Heinkenschloss

  • Feb. 18, 2020

19

slide-22
SLIDE 22

Expansion of Quadratic-Bilinear Systems

◮ Expand around u0(t) = 0 (corresponding state y0(t) = 0). ◮ Assume control u0(t) + au(t) for scalar a > 0. Corresponding state is of form y(t) =

  • ℓ=1

aℓyℓ(t) and output is z(t) = cT y(t) + d

  • u0(t) + au(t)
  • =

  • ℓ=1

aℓcT yℓ(t) + adu(t). ◮ States yℓ in expansion satisfy E d

dty1(t) = Ay1(t) + bu(t),

ℓ = 1, E d

dty2(t) = Ay2(t) + G

  • y1(t), y1(t)
  • + Ny1(t)u(t),

ℓ = 2, E d

dty3(t) = Ay3(t) + G

  • y1(t), y2(t)
  • + G
  • y2(t), y1(t)
  • + Ny2(t)u(t),

ℓ = 3, . . . Given y1, . . . , yℓ−1, ℓ-th equation is linear in yℓ.

Matthias Heinkenschloss

  • Feb. 18, 2020

20

slide-23
SLIDE 23

◮ If E is invertible (assume E = I to shorten notation) solutions are y1(t) = t eAτ1bu(t − τ1)dτ1, y2(t) = t eAτ2[G

  • y1(t − τ2), y1(t − τ2)
  • + Ny1(t − τ2)u(t − τ2)]dτ2,

. . . ◮ Output - here truncated at ℓ = 2: z(t) =

  • ℓ=1

aℓcT yℓ(t) + adu(t) ≈

2

  • ℓ=1

aℓcT yℓ(t) + adu(t), = t h1(τ1)u(t − τ1)dτ1 + t t−τ2 h2(τ1, τ2)u(t − τ1 − τ2)u(t − τ2)dτ1dτ2 + t t−τ3 t−τ3 h3(τ1, τ2, τ3)u(t − τ1 − τ3)u(t − τ2 − τ3)dτ1dτ2dτ3, with kernels h1(τ1) = cT eAτ1b, h2(τ1, τ2) = cT eAτ2NeAτ1b, h3(τ1, τ2, τ3) = cT eAτ3G

  • eAτ2b, eAτ1b
  • .

Matthias Heinkenschloss

  • Feb. 18, 2020

21

slide-24
SLIDE 24

Define Φ(s) = (sE − A)−1. Levels of transfer functions: H0(s) = cT Φ(s)b, H1(s0, s1) = cT Φ(s0)NΦ(s1)b, H2(s0, s1, s2) = cT Φ(s0)G

  • Φ(s1)b, Φ(s2)b
  • ,

. . .

Matthias Heinkenschloss

  • Feb. 18, 2020

22

slide-25
SLIDE 25

Loewner ROM

Frequency data: k = 3¯ k µ(1)

1 , µ(1) 2 , µ(1) 3 , . . . , µ(¯ k) 1 , µ(¯ k) 2 , µ(¯ k) 3 ,

λ(1)

1 , λ(1) 2 , λ(1) 3 , . . . , λ(¯ k) 1 , λ(¯ k) 2 , λ(¯ k) 3 .

Grouped into multi-tuples µ(j) = {(µ(j)

1 ), (µ(j) 1 , µ(j) 2 ), (µ(j) 1 , λ(j) 1 , µ(j) 3 )},

j = 1, . . . , ¯ k, λ(j) = {(λ(j)

1 ), (λ(j) 2 , λ(j) 1 ), (λ(j) 3 , λ(j) 1 , λ(j) 1 )},

j = 1, . . . , ¯ k. Generalized controllability matrix R ∈ Cn×k R =

  • R(1), R(2), · · · , R(¯

k)

∈ Cn×k, R(j) =

  • Φ(λ(j)

1 ) b, Φ(λ(j) 2 ) N Φ(λ(j) 1 ) b, Φ(λ(j) 3 )G

  • Φ(λ(j)

1 )b, Φ(λ(j) 1 )b

  • .

Generalized observability matrix O =

  • O(1)T ,
  • O(2)T , . . .
  • O(¯

k)T T

∈ Ck×n, O(j) =    cT Φ(µ(j)

1 )

cT Φ(µ(j)

1 ) N Φ(µ(j) 2 )

cT Φ(µ(j)

1 ) G

  • Φ(λ(j)

1 )b, Φ(µ(j) 3 )

  .

Matthias Heinkenschloss

  • Feb. 18, 2020

23

slide-26
SLIDE 26

Loewner ROM

◮ Loewner matrix L and shifted Loewner matrix Ls L = −O E R, Ls = −O A R . ◮ Compute (short) SVDs [L Ls] = Y1Σ1X∗

1,

L Ls

  • = Y2Σ2X∗

2,

where Σ1 ∈ Rk×2k, Σ2 ∈ R2k×k, Y1, X2 ∈ Ck×k. ◮ Y, X ∈ Ck×r obtained by selecting first r columns of Y1 and X2. Loewner ROM is

  • E = −Y∗LX = W∗EV,
  • A = −Y∗LsX∗ = W∗EV,
  • G
  • ·, ·
  • = W∗G
  • V·, V ·
  • ,
  • N = W∗NV,
  • b = W∗ b,
  • c = V∗ c.

◮ ROM E, A, . . . , represented by projection, but can compute ROM directly from (sub-)transfer function measurements (no need to explicitly project FOM). ◮ ROM interpolates (sub-)transfer functions.

Matthias Heinkenschloss

  • Feb. 18, 2020

24

slide-27
SLIDE 27

◮ Convergence of expansion for a ∈ (0, amax) proven for finite dimensional, so-called linear analytic state equations (quadratic bilinear systems are a special case) in, e.g., [Rugh, 1981]. ◮ Proof applied to finite dimensional quadratic-bilinear systems would be useful (characterization of amax), as well as extension to PDEs. ◮ Error estimate for truncation needed. Truncation at ℓ = 3 somewhat arbitrary, but seems to work well in practice. Approach can be extended to truncation at ℓ > 3. ◮ Previous derivation assumes E is invertible. Not true for Navier-Stokes. Use (discrete) Leray projections [Heinkenschloss et al., 2008] to extend technique (in progress).

Matthias Heinkenschloss

  • Feb. 18, 2020

25

slide-28
SLIDE 28

Burger’s Equation

Viscosity ν > 0, parameters σ0 ≤ 0, σ1 ≥ 0 (ν = 0.01, σ0 = 0, σ1 = 0.1), ∂ ∂ty(x, t) − ν ∂2 ∂x2 y(x, t) + y(x, t) ∂ ∂xy(x, t) = 0, x ∈ (0, 1), t ∈ (0, T), ν ∂ ∂xy(x, 0) + σ0y(0, t) = u(t), t ∈ (0, T), ν ∂ ∂xy(x, 1) + σ1y(1, t) = 0, t ∈ (0, T), y(x, 0) = 0, x ∈ (0, 1). Input u; output z(t) = 1 y(x, t) dx. Finite element semi-discretization leads to E d

dty(t) = Ay(t) + G(y(t), y(t)) + Ny(t)u(t) + bu(t),

t ∈ (0, T), z(t) = cT y(t), t ∈ (0, T), y(0) = 0.

Matthias Heinkenschloss

  • Feb. 18, 2020

26

slide-29
SLIDE 29

Computations with piecewise linear FEM, FOM size n = 257. Use the input u(t) = 0.1 sin(4πt) in simulations. Frequency data: 300 logarithmically spaced points between 1 and 103 and their complex conjugates (will be able to construct real ROMs). Split into µ1, . . . µk, λ1, . . . λk ∈ C, k = 300. Normalized singular values of the Loewner matrices [L Ls] and L Ls

  • .

ROM size r = 22 is chosen to be smallest r with σr/σ1 > 10−11.

Matthias Heinkenschloss

  • Feb. 18, 2020

27

slide-30
SLIDE 30

Output of FOM and Loewner ROM Error between FOM and Loewner ROM outputs State computed with FOM State computed with Loewner ROM Instabilities may arise because Loewner ROM corresponds to Petrov-Galerkin Projection V = W.

Matthias Heinkenschloss

  • Feb. 18, 2020

28

slide-31
SLIDE 31

◮ Can show existence of constant C > 0 that depends on T, ν > 0, but not on u such that yW (0,T ) + yL∞ ≤ C

  • 1 + uL2(0,T )
  • .

(∗) ◮ Similar estimate holds for any Galerkin approximation. ◮ However, Loewner ROM is a Petrov-Galerkin approximation. Can lead to instabilities - as on previous slide. ◮ Use Loewner ROM projection matrices V, W ∈ Rn×r, construct [V, W] ∈ Rn×2r (orthogonalize), and use this matrix to construct Galerkin projection ROM. This ROM solution obeys estimate like (*).

Matthias Heinkenschloss

  • Feb. 18, 2020

29

slide-32
SLIDE 32

Use Projection Matrix [V, W]

Output of FOM and Loewner ROM Error between FOM and Loewner ROM outputs State computed with FOM State computed with Loewner ROM Excellent agreement between FOM and Loewner ROM (with projection matrix [V, W])

Matthias Heinkenschloss

  • Feb. 18, 2020

30

slide-33
SLIDE 33

Conclusions: ◮ Loewner framework to construct ROM directly from data.

◮ Theoretically, Loewner is a projection based, interpolatory ROM ◮ Can do so for LTI systems. ◮ Can be extended to quadratic-bilinear systems, but not clear how data can be directly obtained.

◮ Loewner can be extended to quadratic-bilinear systems, but not clear how data can be directly obtained. ◮ Error estimates for quadratic-bilinear systems. ◮ Stability of Loewner for quadratic-bilinear systems.

Matthias Heinkenschloss

  • Feb. 18, 2020

31

slide-34
SLIDE 34

Literature I

Ahmad, M. I., Benner, P., Goyal, P., and Heiland, J. (2017). Moment-matching based model reduction for Navier-Stokes type quadratic-bilinear descriptor systems. ZAMM Z. Angew. Math. Mech., 97(10):1252–1267, DOI: 10.1002/zamm.201500262, http://dx.doi.org/10.1002/zamm.201500262. Antoulas, A. C. (2005). Approximation of Large-Scale Dynamical Systems, volume 6 of Advances in Design and Control. Society for Industrial and Applied Mathematics (SIAM), Philadelphia, PA, DOI: 10.1137/1.9780898718713, https://doi.org/10.1137/1.9780898718713. Antoulas, A. C., Beattie, C. A., and Gugercin, S. (2020a). Interpolatory Model Reduction, volume 21 of Computational Science & Engineering. Society for Industrial and Applied Mathematics (SIAM), Philadelphia, PA, DOI: 10.1137/1.9781611976083, https://doi.org/10.1137/1.9781611976083.

Matthias Heinkenschloss

  • Feb. 18, 2020

32

slide-35
SLIDE 35

Literature II

Antoulas, A. C., Gosea, I. V., and Heinkenschloss, M. (2018). On the Loewner framework for model reduction of Burgers’ equation. In King, R., editor, Active Flow and Combustion Control 2018, pages 255–270. Springer-Verlag, Berlin, Heidelberg, New York, DOI: 10.1007/978-3-319-98177-2 16, https://doi.org/10.1007/978-3-319-98177-2_16. Antoulas, A. C., Gosea, I. V., and Heinkenschloss, M. (2020b). Data-driven model reduction for a class of semi-explicit DAEs using the Loewner framework. Technical report, Department of Computational and Applied Mathematics, Rice University. Antoulas, A. C., Lefteriu, S., and Ionita, A. C. (2017). Chapter 8: A tutorial introduction to the Loewner framework for model reduction. In Benner, P., Cohen, A., Ohlberger, M., and Willcox, K., editors, Model Reduction and Approximation: Theory and Algorithms, pages 335–376. SIAM, Philadelphia, DOI: 10.1137/1.9781611974829.ch8, https://doi.org/10.1137/1.9781611974829.ch8.

Matthias Heinkenschloss

  • Feb. 18, 2020

33

slide-36
SLIDE 36

Literature III

Beattie, C. and Gugercin, S. (2017). Chapter 7: Model reduction by rational interpolation. In Benner, P., Cohen, A., Ohlberger, M., and Willcox, K., editors, Model Reduction and Approximation: Theory and Algorithms, Computational Science and Engineering, pages 297–334. SIAM, Philadelphia, DOI: 10.1137/1.9781611974829.ch7, https://doi.org/10.1137/1.9781611974829.ch7. Benner, P. and Breiten, T. (2015). Two-sided projection methods for nonlinear model order reduction. SIAM J. Sci. Comput., 37(2):B239–B260, DOI: 10.1137/14097255X, http://dx.doi.org/10.1137/14097255X. Benner, P. and Breiten, T. (2017). Chapter 6: Model order reduction based on system balancing. In Benner, P., Cohen, A., Ohlberger, M., and Willcox, K., editors, Model Reduction and Approximation: Theory and Algorithms, Computational Science and Engineering, pages 261–295, Philadelphia. SIAM, DOI: 10.1137/1.9781611974829.ch6, https://doi.org/10.1137/1.9781611974829.ch6.

Matthias Heinkenschloss

  • Feb. 18, 2020

34

slide-37
SLIDE 37

Literature IV

Gosea, I. V. and Antoulas, A. C. (2018). Data-driven model order reduction of quadratic-bilinear systems.

  • Numer. Linear Algebra Appl., 25(6):e2200, DOI: 10.1002/nla.2200,

http://dx.doi.org/10.1002/nla.2200. Gosea, I. V., Zhang, Q., and Antoulas, A. C. (2020). Preserving the DAE structure in the Loewner model reduction and identification framework. Advances in Computational Mathematics, 46(1):3, DOI: 10.1007/s10444-020-09752-8, https://doi.org/10.1007/s10444-020-09752-8. Gubisch, M. and Volkwein, S. (2017). Chapter 1: Proper Orthogonal Decomposition for linear-quadratic optimal control. In Benner, P., Cohen, A., Ohlberger, M., and Willcox, K., editors, Model Reduction and Approximation: Theory and Algorithms, Computational Science and Engineering, pages 3–64, Philadelphia. SIAM, DOI: 10.1137/1.9781611974829.ch1, https://doi.org/10.1137/1.9781611974829.ch1.

Matthias Heinkenschloss

  • Feb. 18, 2020

35

slide-38
SLIDE 38

Literature V

Haasdonk, B. (2017). Chapter 2: Reduced basis methods for parametrized PDEs - a tutorial introduction for stationary and instationary problems. In Benner, P., Cohen, A., Ohlberger, M., and Willcox, K., editors, Model Reduction and Approximation: Theory and Algorithms, Computational Science and Engineering, pages 65–136. SIAM, Philadelphia, DOI: 10.1137/1.9781611974829.ch2, https://doi.org/10.1137/1.9781611974829.ch2. Heinkenschloss, M., Sorensen, D. C., and Sun, K. (2008). Balanced truncation model reduction for a class of descriptor systems with application to the Oseen equations. SIAM Journal on Scientific Computing, 30(2):1038–1063, DOI: 10.1137/070681910, http://doi.org/10.1137/070681910. Hesthaven, J. S., Rozza, G., and Stamm, B. (2015). Certified Reduced Basis Methods for Parametrized Partial Differential Equations. Springer Briefs in Mathematics. Springer, New York, DOI: 10.1007/978-3-319-22470-1, http://doi.org/10.1007/978-3-319-22470-1.

Matthias Heinkenschloss

  • Feb. 18, 2020

36

slide-39
SLIDE 39

Literature VI

Quarteroni, A., Manzoni, A., and Negri, F. (2016). Reduced Basis Methods for Partial Differential Equations. An Introduction, volume 92 of Unitext. Springer, Cham, DOI: 10.1007/978-3-319-15431-2, https://doi.org/10.1007/978-3-319-15431-2. Rowley, C. W. and Dawson, S. T. M. (2017). Model reduction for flow analysis and control. Annual Review of Fluid Mechanics, 49(1):387–417, DOI: 10.1146/annurev-fluid-010816-060042, https://doi.org/10.1146/annurev-fluid-010816-060042. Rozza, G., Huynh, D. B. P., and Patera, A. T. (2008). Reduced basis approximation and a posteriori error estimation for affinely parametrized elliptic coercive partial differential equations: application to transport and continuum mechanics.

  • Arch. Comput. Methods Eng., 15(3):229–275, DOI:

10.1007/s11831-008-9019-9, http://dx.doi.org/10.1007/s11831-008-9019-9.

Matthias Heinkenschloss

  • Feb. 18, 2020

37

slide-40
SLIDE 40

Literature VII

Rugh, W. J. (1981). Nonlinear System Theory. The Volterra/Wiener approach. Johns Hopkins Series in Information Sciences and Systems. Johns Hopkins University Press, Baltimore, Md. PDF version with corrections (2002) https://sites.google.com/site/wilsonjrugh (accessed Feb 13, 2018).

Matthias Heinkenschloss

  • Feb. 18, 2020

38