Singular perturbations, algorithmic algebra and scalings Sebastian - - PowerPoint PPT Presentation

Singular perturbations, algorithmic algebra and scalings

Sebastian Walcher, RWTH Aachen Symbiont Meeting, Bonn, July 2018

Based on work by and work with: Lena N¨

• then

Alexandra Goeke Christian Lax Eva Zerz

Overview

Objective: Dimension reduction for polynomial and rational ODEs. Motivated by chemical reaction networks, with slow and fast reactions, slow and fast variables (QSS). Mathematics: Singular perturbations. Contents:

◮ Method of choice (and restriction): Classical reduction

theorems by Tikhonov and Fenichel; coordinate-free. (Reliable convergence results, verifiable conditions.)

◮ Tikhonov-Fenichel and QSS. ◮ Tikhonov-Fenichel and scaling approaches.

Singular perturbation reduction

Tikhonov and Fenichel (reminder)

System with small parameter ε in standard form ˙ y1 = εf (y1, y2) + ε2 (. . . ), y1 ∈ D ⊆ Rs, ˙ y2 = g(y1, y2) + ε (. . . ), y2 ∈ G ⊆ Rr Slow time τ = εt : y′

1 = f (y1, y2)+· · · ,

εy′

2 = g(y1, y2)+ε · · · .

Assumptions: (i) Critical manifold

• Z :=
• (y1, y2)T ∈ D × G; g(y1, y2) = 0
• = ∅;

(ii) there exists ν > 0 such that every eigenvalue of D2g(y1, y2), (y1, y2) ∈ Z has real part ≤ −ν.

• Theorem. There exist T > 0 and a neighborhood of

Z in which all solutions converge uniformly to solutions of y′

1 = f (y1, y2),

g(y1, y2) = 0

• n [t0, T] ; t0 > 0 arbitrary.

Ordinary differential equations for CRN

We consider: Parameter dependent ordinary differential equations ˙ x = h(x, π), x ∈ Rn, π ∈ Rm with polynomial (or rational) right hand side. Motivation: Chemical reaction networks with mass action kinetics and constant thermodynamical parameters; spatially homogeneous. Parameters: Rate constants, initial concentrations. Question: Where (and how) do singular perturbations enter the picture?

Tikhonov-Fenichel, coordinate-free: First step

We have ˙ x = h(x, π) versus ˙ y1 = εf (y1, y2) + ε2 (. . . ) ˙ y2 = g(y1, y2) + ε (. . . ) Fact and problems: No a priori separation in slow and fast

• variables. No small parameter ε.

Preliminary step: For suitable (?!) π consider system ˙ x = h(x, π + ερ + · · · ) =: h(0)(x) + ε h(1)(x) + ε2 · · · . Suitability of π: Scenario must be singular; i.e. the vanishing set V(h(0)) contains a submanifold Z with dimension s > 0.

Tikhonov-Fenichel, coordinate-free: Identification

• Proposition. Assume dim Z = s > 0. Then

˙ x = h(0)(x) + ε h(1)(x) + ε2 . . . admits a coordinate transformation into standard form and subsequent Tikhonov-Fenichel reduction near every point of Z if and only if (i) rank D h(0)(x) = r := n − s for all x ∈ Z; (ii) for each x ∈ Z there exists a direct sum decomposition Rn = Ker D h(0)(x) ⊕ Im D h(0)(x); (iii) for each x ∈ Z the nonzero eigenvalues of D h(0)(x) have real parts ≤ −ν < 0. Remaining problem: Direct application (via explicitly computable coordinate transformation) may be impossible.

Tikhonov-Fenichel, coordinate-free: Reduction

x′ = ε−1 h(0)(x) + h(1)(x) + . . . with a ∈ Z ⊆ V( h(0)), satisfying conditions (i), (ii) und (iii). Decomposition: There is a Zariski-open neighborhood Ua of a such that

• h(0)(x) = P(x)µ(x),

with µ(x) ∈ R(x)r×1, P(x) ∈ R(x)n×r, rank P(a) = r, rank Dµ(a) = r, and (wlog) V( h(0)) ∩ Ua = V(µ) ∩ Ua = Z. Reduction: The system x′ =

• In − P(x)A(x)−1Dµ(x)
• h(1)(x),

with A(x) := Dµ(x)P(x) is defined on Ua and admits Z as invariant set. The restriction to Z yields the reduction from Tikhonov’s theorem as ε → 0.

Tikhonov-Fenichel, coordinate-free: Parameters

Definition: We call π a Tikhonov-Fenichel parameter value (TFPV) for dimension s (1 ≤ s ≤ n − 1) of ˙ x = h(x, π) if the following hold: (i) The vanishing set V(h(·, π)) of x → h(x , π) contains a component Y of dimension s. (ii) There is x0 ∈ Y and neighborhood Z of x0 in Y such that rank D1h(x, π) = n − s and Rn = Ker D1h(x, π) ⊕ Im D1h(x, π), for all x ∈ Z. (iii) The nonzero eigenvalues of D1h(x0, π) have real parts < 0. Note: Conditions simply were copied from above. Therefore reduction works for small perturbations π + ερ + · · · .

Determining TFPV (algorithmic starting point)

• Proposition. If

π is a TFPV for dimension s (1 ≤ s ≤ n − 1) of ˙ x = h(x, π) then there exists x0 ∈ Rn such that: (i) h(x0, π) = 0; (ii) for each k > r := n − s all k × k minors of the Jacobian D1h(x0, π) vanish. Observation: More equations for x than variables. Enter algorithmic algebra (elimination theory). More theory:

• Theorem. The TFPV of a polynomial (or rational) system

˙ x = h(x, π) with nonnegative parameters form a semi-algebraic subset of Rm.

Standard example: Michaelis-Menten (irreversible)

Reaction network E + S

k1

⇋

k−1 C k2

⇀ E + P Differential equation system for concentrations: ˙ s = −k1es + k−1c, ˙ c = k1es − (k−1 + k2)c, ˙ e = −k1es + (k−1 + k2)c; typical initial values s(0) = s0, c(0) = 0, e(0) = e0. With stoichiometry (linear first integral e + c): ˙ s = − k1e0s + (k1s + k−1)c, ˙ c = k1e0s − (k1s + k−1 + k2)c.

Michaelis-Menten: Determine TFPV

System ˙ s = − k1e0s + (k1s + k−1)c, ˙ c = k1e0s − (k1s + k−1 + k2)c with Jacobian determinant d = k1k2(e0 − c). Eliminate s and c from the three equations (manageable by hand). Result: A TFPV ( e0, k1, k−1, k2) satisfies

• e0

k2 k1 = 0. Small perturbations yield (all) relevant cases:

    

εe∗

• k1
• k−1
• k2

     or     

• e0

εk∗

1

• k−1
• k2

     or     

• e0
• k1
• k−1

εk∗

2

     or     

• e0
• k1

εk∗

−1

εk∗

2

    

Michaelis-Menten; some reductions

◮ Michaelis-Menten with little enzyme e0 = εe∗ 0: Familiar result. ◮ Michaelis-Menten with slow product formation:

˙ s = − k1e0s + (k1s + k−1)c ˙ c = k1e0s − (k1s + k−1)c − εk∗

2c.

Decomposition h(0) = P · µ with P =

• 1

−1

• , µ = k1e0s − (k1s + k−1)c.

Reduced equation (on Z = V(µ)):

• s′

c′

• =

1 k1(e0 − c) + k1s + k−1

• ∗

k1s + k−1 ∗ k1(e0 − c)

• ·
• −k∗

2c

• .

Further example: TFPV for competitive inhibition

Michaelis-Menten network with inhibitor (binding reversibly to enzyme) leads to three dimensional system ˙ s = k−1c1 − k1s(e0 − c1 − c2), ˙ c1 = k1s(e0 − c1 − c2) − (k−1 + k2)c1, ˙ c2 = k3(e0 − c1 − c2)(i0 − c2) − k−3c2 with Jacobian determinant d(x, π) = −k1k2(e0 − c1 + c2)(k−3 + k3(i0 + e0) − k3(2c2 − c1)). Elimination ideal has radical I = e0k1k2k−3(k2

3(e0 − i0)2 + k−3(k−3 + 2k3(e0 + i0))

with single generator.

Tihkonov-Fenichel approach: Review

◮ The approach is (in principle) universal. (Restriction to two

time scales due to theoretical foundation.)

◮ Defining equations and inequalities algorithmically accessible

(in weak sense: termination guaranteed but no criterion).

◮ Methods in current use not quite adequate (algebra rather

than real algebra).

◮ Obvious practical problem: Feasibility.

More restrictive searches for slow manifolds (e.g. by scaling) remain relevant, from scientific (chemistry) and from pragmatic perspective.

Quasi-steady reduction

“Classical” QSS reduction

Quasi-steady state (QSS) reduction of a parameter dependent system ˙ x = h(x, π).

◮ Loose description: Designate some variables as being in

“quasi-equilibrium” (zero rate of change) after short initial

• phase. This yields a system of algebraic equations which are

used to eliminate these variables. (Heuristical procedure.)

◮ First explicit statement: Briggs and Haldane (1925) for

Michaelis-Menten: QSS for c.

◮ Obvious: Consistency and validity of this procedure depend

• n conditions for parameters (Briggs and Haldane require

small enzyme concentration e0).

Tikhonov-Fenichel vs.“classical” QSS reduction

◮ Reductions agree for Michaelis-Menten with small enzyme

• concentration. (Heineken, Tsuchiya and Aris 1967: Possibly

the first mathematically rigorous application of Tikhonov to a reaction network.)

◮ Concepts are frequently conflated, but should not be. ◮ QSS reduction may yield incorrect results. ◮ Singular perturbation reductions do not only target slow and

fast species: For instance, slow and fast reactions are equally relevant.

An incorrect QSS reduction

Example: Irreversible Michaelis-Menten equation ˙ s = − k1e0s + (k1s + k−1)c ˙ c = k1e0s − (k1s + k−1)c − εk∗

2c

with slow product formation; k2 = εk∗

2.

Tikhonov-Fenichel reduction on critical manifold Z (given by k1e0s − (k1s + k−1)c = 0): ˙ s = − k2k1e0s (k1s + k−1) k−1e0 + (k1s + k−1)2 . QSS reduction for complex (lowest order in ε): ˙ s = − k2k1e0s k1s + k−1 + k2 = − k2k1e0s k1s + k−1 + · · · These differ significantly (and QSS is wrong)!

QSS and Tikhonov: A positive result

• Proposition. Let

π be a TFPV of the parameter dependent system ˙ x = h(x, π), with critical manifold a coordinate subspace, given by xs+1 = γs+1, . . . , xn = γn (up to relabelling). Then the Tikhonov-Fenichel reduction of ˙ x = h(x, π + ερ + . . .) = h(x, π) + εD2h(x, π)ρ + . . . and the classical quasi-steady state reduction of this equation with respect to xs+1, . . . , xn agree up to first order in ε. Consequence: Tikhonov ensures correctness of the QSS reduction. Can one do better than this?

Singular perturbations and scaling

Scaling: Simplest version

Given system with “small parameter” ε: ˙ x = h(x, ε) = h(0)(x) + ε h(1)(x) + · · · . For a partitioning x =

• z

y

• with z ∈ Rr, y ∈ Rs, r > 1, s > 1, r + s = n,

rewrite system with initial value x0 = (z0, y0)tr as ˙ z = f (z, y, ε) = f0(z, y) + εf1(z, y) + ε2 · · · , ˙ y = g(z, y, ε) = g0(z, y) + εg1(z, y) + ε2 · · · . (1) An asymptotically degenerate scaling is a linear transformation

• z

y

• =
• z

εy∗

• .

Scaling: A historical example

Heineken, Tsuchiya und Aris: Reduction of ˙ s = −k1e0s + (k1s + k−1)c, ˙ c = k1e0s − (k1s + k−1 + k2)c. Small enzyme concentration, translated to e0 = εe∗

0, ε → 0.

Degenerate scaling : Set c = εc∗ to obtain system ˙ s = ε(−k1se∗

0 + (k1s + k−1)c∗),

˙ c∗ = k1s − (k1s + k−1 + k2)c∗ in standard form. Applying Tikhonov yields familiar reduced equation ˙ s = −k1e0s/(k1s + k−1 + k2). Warning notice: This procedure works due to special (consistency) properties of the system.

Scaling: Consistency

Replacing y by εy∗ in (1) and rearranging yields ˙ z = f0(z, 0) + ε(f1(z, 0) + D2f0(z, 0)y∗) + ε2 · · · ˙ y∗ = ε−1g0(z, 0) + (g1(z, 0) + D2g0(z, 0)y∗) + ε · · · (2) with initial values z(0) = z0, y∗(0) = y∗

0 := ε−1y0.

(a) We call the scaling transformation y = εy∗ (and system (2)) locally consistent if g0(z, 0) = 0 for all z ∈ U. We call the scaling transformation locally Tikhonov consistent if it is locally consistent and furthermore f0(z, 0) = 0 for all z ∈ U (cue non-isolated stationary points). (b) We call the scaling y = εy∗ initial value consistent if the initial value for y satisfies y0 = εy∗

0 , with fixed y∗ 0 .

• Remark. If these consistency conditions do not hold then “singular

perturbation reductions” may be meaningless or incorrect.

Initial value consistency: Example

Two-dimensional linear system ˙ z = −εaz + by ˙ y = −cy with constants a ≥ 0, b, and c > 0, and initial values z0 resp. y0. Scale y = εy∗ and pass to slow time τ: z′ = −az + by∗ (y∗)′ = −ε−1cy∗ with initial values z(0) = z0 and y∗(0) = ε−1y0. The scaled system seems amenable to Tikhonov, but “escaping” initial value y∗(0) precludes this. (Verification straightforward.)

Scalings and reduction: Standard case

◮ Whenever the scaling y = εy∗ satisfies local Tikhonov

consistency and initial value consistency then system (2) locally has the form ˙ z = + ε(f1(z, 0) + F0(z, 0)y∗) + · · · ˙ y∗ = G0(z, 0)y∗ + g1(z, 0) + ε · · · with suitable matrix-valued functions F0 and G0 and initial conditions z(0) = z0, y∗(0) = y∗

0 . ◮ Standard case (eigenvalue conditions and invertibility for G0)

yields reduced system (agreeing with QSS reduction) z′ = f1(z, 0) − F0(z, 0)G0(z, 0)−1g1(z, 0), z(0) = z0. Note: Direct use of coordinate-free Tikhonov-Fenichel is easier. But there is a nonstandard case (G0 not invertible).

Scalings and reduction: Nonstandard example

Michaelis-Menten with slow degradation of free enzyme: ˙ s = −k1es + k−1c ˙ e = −k1es + (k−1 + k2)c − εδ∗e ˙ c = k1es − (k−1 + k2)c Scale e = εe∗, c = εc∗; then G0 =

• −k1s

k−1 + k2 k1s −(k−1 + k2)

• is not invertible. Reduced system

s′ = −k1e∗s + k−1c∗ e∗′ = −k1e∗

d

• −k1e∗s + k−1c∗ + δ∗(k−1+k2)

k1

• c∗′

=

k1e∗ d

(−k1e∗s + k−1c∗ − δ∗s)

• n 2D manifold defined by µ := (k1e∗s − (k−1 + k2)c∗) = 0.

Determining consistent scalings: Example

Michaelis Menten (again!), no slow reactions: ˙ s = −k1es + k−1c ˙ e = −k1es + (k−1 + k2)c ˙ c = k1es − (k−1 + k2)c

◮ Local Tikhonov consistency: Find partitioning

x =

• z

y

• with

h(0)(z, 0) = 0; y = εy∗.

◮ One must have c = εc∗. ◮ Then two choices e = εe∗ (usual reduction) or s = εs∗ (trivial

reduced equation).

◮ Note: Initial value consistency enters e.g. for c(0) and e(0).

Scalings: Review and outlook

◮ Simple scaling approach not necessary for the purpose of

Heineken et al.

◮ But indispensable for nonstandard cases. Class of

applications: Reaction-transport systems, e.g. M-M with diffusion.

◮ Consistency conditions may be employed to find scalings. ◮ Consistency conditions (in particular initial value consistency)

should provide further conditions for tropical geometry approach.

Thank you for your attention!