{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"slideshow": {
"slide_type": "skip"
}
},
"outputs": [],
"source": [
"from __future__ import division\n",
"from warnings import filterwarnings\n",
"filterwarnings('ignore')\n",
"\n",
"from matplotlib import pyplot as plt\n",
"import numpy as np\n",
"\n",
"%matplotlib inline"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"## AERSP 597 - Machine Learning in Aerosapce Engineering\n",
"## Lecture 25, Dynamical System: Intro to Koopman Theory\n",
"### Instructor: Daning Huang"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "-"
}
},
"source": [
"$$\n",
"\\newcommand{\\bR}{\\mathbb{R}}\n",
"\\newcommand{\\vb}{\\mathbf{b}}\n",
"\\newcommand{\\vf}{\\mathbf{f}}\n",
"\\newcommand{\\vg}{\\mathbf{g}}\n",
"\\newcommand{\\vu}{\\mathbf{u}}\n",
"\\newcommand{\\vv}{\\mathbf{v}}\n",
"\\newcommand{\\vw}{\\mathbf{w}}\n",
"\\newcommand{\\vx}{\\mathbf{x}}\n",
"\\newcommand{\\vy}{\\mathbf{y}}\n",
"\\newcommand{\\vA}{\\mathbf{A}}\n",
"\\newcommand{\\vB}{\\mathbf{B}}\n",
"\\newcommand{\\vG}{\\mathbf{G}}\n",
"\\newcommand{\\vU}{\\mathbf{U}}\n",
"\\newcommand{\\vV}{\\mathbf{V}}\n",
"\\newcommand{\\vW}{\\mathbf{W}}\n",
"\\newcommand{\\vX}{\\mathbf{X}}\n",
"\\newcommand{\\vtt}{\\boldsymbol{\\theta}}\n",
"\\newcommand{\\vtf}{\\boldsymbol{\\phi}}\n",
"\\newcommand{\\vty}{\\boldsymbol{\\psi}}\n",
"\\newcommand{\\vtx}{\\boldsymbol{\\xi}}\n",
"\\newcommand{\\vtF}{\\boldsymbol{\\Phi}}\n",
"\\newcommand{\\vtG}{\\boldsymbol{\\Gamma}}\n",
"\\newcommand{\\vtS}{\\boldsymbol{\\Sigma}}\n",
"\\newcommand{\\vtL}{\\boldsymbol{\\Lambda}}\n",
"\\newcommand{\\vtX}{\\boldsymbol{\\Xi}}\n",
"\\newcommand{\\vtY}{\\boldsymbol{\\Psi}}\n",
"\\newcommand{\\cB}{\\mathcal{B}}\n",
"\\newcommand{\\cF}{\\mathcal{F}}\n",
"\\newcommand{\\cK}{\\mathcal{K}}\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"## TODAY: Dynamical System - V\n",
"\n",
"* Koopman operator\n",
"* Extended DMD\n",
"* Control application\n",
"\n",
"### References:\n",
"\n",
"* [DMSC] Part 3\n",
"* [Brunton2022] Modern Koopman Theory for Dynamical Systems, 2022"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"In L23, we have presented three underpinning principles for data-driven modeling of dynamical systems. The last one is that\n",
"\n",
"> **Global linearization**: A **finite**-dimensional nonlinear system can be represented by an **infinite**-dimensional linear system"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"We have shown a special example where the linear system is actually finite-dimensional,\n",
"$$\n",
"\\left\\{\\begin{array}\n",
"\\ \\dot x_1 = \\mu x_1 \\\\\n",
"\\dot x_2 = \\lambda (x_2-x_1^2)\n",
"\\end{array}\\right.\n",
"\\Rightarrow\n",
"\\left\\{\\begin{array}\n",
"\\ \\dot y_1 = \\mu y_1 \\\\\n",
"\\dot y_2 = \\lambda (y_2-y_3) \\\\\n",
"\\dot y_3 = 2\\mu y_3\n",
"\\end{array}\\right.\n",
"$$\n",
"where\n",
"$$\n",
"y_1=x_1,\\quad y_2=x_2,\\quad y_3=x_1^2\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"Here is one more example of PDE - 1D viscous Burgers' equation (model problem for shock waves)\n",
"$$\n",
"\\frac{\\partial u}{\\partial t}+u\\frac{\\partial u}{\\partial x} = \\nu \\frac{\\partial^2 u}{\\partial x^2},\\quad u(x,0)=f(x)\n",
"$$\n",
"\n",
"\n",
"\n",
"From Wikipedia, [Cole-Hopf Transformation](https://en.wikipedia.org/wiki/Burgers%27_equation#Viscous_Burgers'_equation)"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"Choosing\n",
"$$\n",
"u = -2\\nu\\frac{1}{\\phi}\\frac{\\partial \\phi}{\\partial x}\n",
"$$\n",
"one gets a linear PDE that can be solved easily\n",
"$$\n",
"\\frac{\\partial \\phi}{\\partial t} = \\nu\\frac{\\partial^2 \\phi}{\\partial x^2}\n",
"$$\n",
"(given the correct boundary conditions)"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"### However...\n",
"\n",
"In general, the \"global linearization\" in finite dimensions does not always work.\n",
"\n",
"For example, if we modify the first example as\n",
"$$\n",
"\\begin{align}\n",
"\\dot x_1 &= \\mu x_1 \\\\\n",
"\\dot x_2 &= \\lambda \\color{red}{(x_2^2-x_1)}\n",
"\\end{align}\n",
"$$\n",
"Then the previous strategy will result in a cascading problem, i.e. infinite $y_i$'s will be needed to make the transformed system linear,\n",
"$$\n",
"\\begin{align}\n",
"\\dot y_1 &= \\mu y_1 \\\\\n",
"\\dot y_2 &= \\lambda (y_3-y_1) \\\\\n",
"\\dot y_3 &= 2 \\lambda (y_4-y_5) \\\\\n",
"\\dot y_4 &= \\cdots \\\\\n",
"\\vdots &= \\cdots\n",
"\\end{align}\n",
"$$\n",
"where\n",
"$$\n",
"y_1=x_1,\\quad y_2=x_2,\\quad y_3=x_2^2,\\quad y_4=x_2^3,\\quad y_5=x_1x_2,\\quad \\cdots\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"## Koopman Operator (Discrete time)\n",
"\n",
"Now let's formalize the \"global linearization\" procedure. Let's first start with the discrete time case,\n",
"$$\n",
"\\begin{align}\n",
"\\mbox{System dynamics: }&\\ \\vx_{k+1}=\\vf(\\vx_k) \\\\\n",
"\\mbox{Observation: }&\\ \\vg_k=\\vg(\\vx_k) \\\\\n",
"\\mbox{Initial condition: }&\\ \\vg_0=\\vg(\\vx_0)\n",
"\\end{align}\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"Next, let's first just consider one scalar observation $g(\\vx)$ -\n",
"\n",
"A Koopman operator $\\cK_t$ is infinite-dimensional, linear, and,\n",
"$$\n",
"\\cK_t g = g\\circ\\vf\n",
"$$\n",
"\n",
"> In the notation above, $\\cK_t$ is a \"functional\" that takes in a function ($g$) and returns another function. $\\circ$ means function composition.\n",
"\n",
"Essentially, $\\cK_t$ does the following, observable **linearizing**:\n",
"$$\n",
"\\cK_t g(\\vx_k) = g(\\vf(\\vx_k)) = g(\\vx_{k+1})\n",
"$$\n",
"\n",
"> The first step is by definition of $\\cK_t$ and the second step is due to system dynamics"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"\n",
"\n",
"Figure Credit: Gueho, D., et al., DOI: 10.2514/6.2022-0655"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"$\\cK_t$ is linear, so\n",
"$$\n",
"\\cK_t(c_1 g_1+c_2 g_2) = c_1 \\cK_t g_1+c_2 \\cK_t g_2\n",
"$$\n",
"for any scalars $c_1,c_2$ and any functions $g_1,g_2$.\n",
"\n",
"And for linear transformations (even infinite dimensional), one can always define an eigenvalue problem\n",
"$$\n",
"\\cK_t\\phi_i = \\tilde{\\lambda}_i\\phi_i,\\quad i=0,1,\\cdots\n",
"$$\n",
"where $\\tilde{\\lambda}_i$ are Koopman eigenvalues, and $\\phi_i$ are Koopman **eigenfunctions**.\n",
"\n",
"> There are infinite number of eigenfunctions, since $\\cK_t$ is infinite-dim"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"From a perspective of discrete-time system, one should expect\n",
"$$\n",
"\\tilde{\\lambda}_i\\phi_i(\\vx_k) = \\cK_t\\phi_i(\\vx_k) = \\phi(\\vf(\\vx_k)) = \\phi(\\vx_{k+1})\n",
"$$\n",
"\n",
"But the above \"trivial\" relation leads to an interesting **lattice** property - consider the product of two eigenfunctions $\\phi(\\vx)=\\phi_i(\\vx)\\phi_j(\\vx)$\n",
"$$\n",
"\\begin{align}\n",
"\\cK_t(\\phi(\\vx)) &= \\cK_t(\\phi_i(\\vx)\\phi_j(\\vx)) \\\\\n",
"&= \\phi_i(\\vf(\\vx))\\phi_j(\\vf(\\vx)) \\\\\n",
"&= \\tilde{\\lambda}_i\\tilde{\\lambda}_j\\phi_i(\\vx)\\phi_j(\\vx)\n",
"\\end{align}\n",
"$$\n",
"\n",
"> The product of Koopman eigenfunctions is also a Koopman eigenfunction! Same for eigenvalues"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"## Koopman Mode Decomposition\n",
"Just like eigenvectors in finite-dim space, eigenfunctions in infinite-dim space form a basis for representing other functions.\n",
"\n",
"For a scalar observation function $g_i(\\vx)$\n",
"$$\n",
"g_i(\\vx) = \\sum_{j=0}^\\infty \\xi_{ij}\\phi_j(\\vx) =\n",
"\\begin{bmatrix} \\xi_{11} & \\xi_{12} & \\xi_{13} & \\cdots \\end{bmatrix}\n",
"\\begin{bmatrix} \\phi_1(\\vx) \\\\ \\phi_2(\\vx) \\\\ \\phi_3(\\vx) \\\\ \\vdots \\end{bmatrix}\n",
"$$\n",
"where $\\xi_{ij}$ are scalar coefficients."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"When there are multiple observables, each observable can be expanded likewise,\n",
"$$\n",
"\\begin{bmatrix} g_1(\\vx) \\\\ g_2(\\vx) \\\\ \\vdots \\\\ g_m(\\vx) \\end{bmatrix} =\n",
"\\begin{bmatrix}\n",
"\\xi_{11} & \\xi_{12} & \\xi_{13} & \\cdots \\\\\n",
"\\xi_{21} & \\xi_{22} & \\xi_{23} & \\cdots \\\\\n",
"\\vdots & \\vdots & \\vdots & \\cdots \\\\\n",
"\\xi_{m1} & \\xi_{m2} & \\xi_{m3} & \\cdots\n",
"\\end{bmatrix}\n",
"\\begin{bmatrix} \\phi_1(\\vx) \\\\ \\phi_2(\\vx) \\\\ \\phi_3(\\vx) \\\\ \\vdots \\end{bmatrix}\n",
"$$\n",
"Or in short,\n",
"$$\n",
"\\vg(\\vx) = \\sum_{i=0}^\\infty \\phi_i(\\vx)\\vtx_i\n",
"$$\n",
"where $\\vtx_i$ are columns of the coefficient matrix and called the **Koopman modes**.\n",
"\n",
"> This is the **Koopman Mode Decomposition**"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"Then, the time response is represented as\n",
"$$\n",
"\\vg(\\vx_{k+1}) = \\cK_t\\vg(\\vx_k) = \\sum_{i=0}^\\infty \\cK_t\\phi_i(\\vx_k)\\vtx_i = \\sum_{i=0}^\\infty \\tilde{\\lambda}_i\\phi_i(\\vx_k)\\vtx_i\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"One can further go backwards in time\n",
"$$\n",
"\\begin{align}\n",
"\\vg(\\vx_{k+1}) &= \\cK_t^2\\vg(\\vx_{k-1}) = \\sum_{i=0}^\\infty \\cK_t^2\\phi_i(\\vx_{k-1})\\vtx_i = \\sum_{i=0}^\\infty \\tilde{\\lambda}_i^2\\phi_i(\\vx_{k-1})\\vtx_i \\\\\n",
"&= \\vdots \\\\\n",
"&= \\sum_{i=0}^\\infty \\tilde{\\lambda}_i^{k+1}\\phi_i(\\vx_0)\\vtx_i\n",
"\\end{align}\n",
"$$\n",
"where the initial condition is expanded as\n",
"$$\n",
"\\vg(\\vx_0) = \\sum_{i=0}^\\infty \\phi_i(\\vx_0)\\vtx_i\n",
"$$\n",
"\n",
"> Looking familiar?"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"#### Recall\n",
"If we collect the observables as a time sequence, and apply DMD to identify a discrete linear system, $$\n",
"\\vg_{k+1} = \\vA \\vg_k\n",
"$$\n",
"where $\\vA\\in\\bR^{m\\times m}$ has eigenvalues $\\lambda_i$ and eigenvectors $\\vtf_i$ (i.e. the DMD modes).\n",
"\n",
"Then the solution can be written as,\n",
"$$\n",
"\\vg_{k+1}=\\sum_{i=1}^m b_i \\lambda_i^{k+1} \\vtf_i\n",
"$$\n",
"where the initial condition is expanded as\n",
"$$\n",
"\\vg_0 = \\sum_{i=1}^m b_i \\vtf_i\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"It is clear there is a one-to-one correspondence\n",
"+ DMD modes $\\vtf_i$ v.s Koopman modes $\\vtx_i$\n",
"+ DMD eigenvalues $\\lambda_i$ v.s. Koopman eigenvalues $\\tilde{\\lambda}_i$\n",
"+ DMD modal coordinates $b_i$ v.s. Koopman eigenfunctions $\\phi_i(\\vx_0)$ (or more generally for any $\\vx_k$)\n",
"\n",
"Therefore DMD is an approximation of the Koopman operator via truncation.\n",
"> This explains why DMD works for many nonlinear problems (having a low-dim linear structure, such as a limit cycle)."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"## Extended DMD\n",
"\n",
"The connection between DMD and Koopman also provides a means for numerical implementation of Koopman mode decomposition."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"First, construct the data matrices as\n",
"$$\n",
"\\vG=[\\vg_0,\\vg_1,\\cdots,\\vg_N],\\quad \\vG^1=[\\vg_1,\\vg_2,\\cdots,\\vg_{N+1}]\n",
"$$\n",
"\n",
"From standard DMD, we can find the system matrix, together with its eigendecomposition,\n",
"$$\n",
"\\vA = \\vG^1\\vG^+ \\approx \\vtF\\vtL\\vtY^H\n",
"$$\n",
"where $+$ is pseudo-inverse.\n",
"\n",
"> Based on previous discussion, $\\vtL$ are simultaneously the DMD and Koopman eigenvalues"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"Next, let's recover the Koopman eigenfunctions.\n",
"\n",
"An eigenfunction and its eigenvalue should satisfy\n",
"$$\n",
"\\begin{bmatrix} \\lambda\\phi(\\vx_0) \\\\ \\lambda\\phi(\\vx_1) \\\\ \\vdots \\\\ \\lambda\\phi(\\vx_N) \\end{bmatrix}\n",
"=\\begin{bmatrix} \\phi(\\vx_1) \\\\ \\phi(\\vx_2) \\\\ \\vdots \\\\ \\phi(\\vx_{N+1}) \\end{bmatrix}\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"Think the observables as basis functions that approximate an eigenfunction,\n",
"$$\n",
"\\phi(\\vx) \\approx \\sum_{i=1}^m g_i(\\vx)u_i = \\vg(\\vx)^T\\vu\n",
"$$\n",
"Then\n",
"$$\n",
"\\begin{bmatrix} \\phi(\\vx_1) \\\\ \\phi(\\vx_2) \\\\ \\vdots \\\\ \\phi(\\vx_{N+1}) \\end{bmatrix} =\n",
"\\begin{bmatrix} \\vg(\\vx_1)^T \\\\ \\vg(\\vx_2)^T \\\\ \\vdots \\\\ \\vg(\\vx_{N+1})^T \\end{bmatrix}\\vu = (\\vG^1)^T\\vu\n",
"$$\n",
"and doing the same for the left-hand side\n",
"$$\n",
"\\lambda\\vG^T\\vu = (\\vG^1)^T\\vu\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"Thinking again in least squares, we get an eigenvalue problem\n",
"$$\n",
"\\lambda\\vu = (\\vG^T)^+(\\vG^1)^T\\vu = \\vA^T\\vu\n",
"$$\n",
"So $\\vu$ are the left eigenvectors of $\\vA$, i.e. $\\vtY$. In other words, EDMD approximates the Koopman eigenfunctions as $\\vg(\\vx)^T\\vtY$."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"#### Side note\n",
"\n",
"Why can we truncate \"safely\" an infinite-dimensional system? There are several arguments\n",
"+ From functional analysis, if one approximates a function by, say eigenfunctions $\\{\\phi_i\\}_{i=1}^\\infty$, the coefficients always decay as $i\\rightarrow\\infty$.\n",
"+ From dynamical systems, the trajectory is more likely residing on a bounded low-dimensional structure (limit cycle, strange attractor, etc.), so a handful eigenfunctions might be sufficient to characterize the structure.\n",
"+ The lattice property of $\\cK_t$ indicates that, in at least some cases, the infinite eigenpairs might be simply generated by a finite set of eigenpairs."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"## Koopman Operator (Continuous time)\n",
"\n",
"Now consider a dynamical system that is sufficiently smooth,\n",
"$$\n",
"\\begin{align}\n",
"\\mbox{System dynamics: }&\\ \\dot{\\vx}=\\vf(\\vx) \\\\\n",
"\\mbox{Observation: }&\\ \\vg=\\vg(\\vx) \\\\\n",
"\\mbox{Initial condition: }&\\ \\vx(0)=\\vx_0\n",
"\\end{align}\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"One can seek the **Koopman generator**\n",
"$$\n",
"\\dot{g}=\\cK g\n",
"$$\n",
"where note that we drop the subscript $t$.\n",
"\n",
"$\\cK$ is essentially the infinitesimal limit of the Koopman operator, so that\n",
"$$\n",
"\\cK g = \\lim_{t\\rightarrow 0} \\frac{\\cK_t g-g}{t} = \\lim_{t\\rightarrow 0} \\frac{g\\circ\\vf-g}{t}\n",
"$$\n",
"\n",
"> In some sense, $\\cK$ and $\\cK_t$ are analogous to $\\vA$ and $\\exp(\\vA t)$ for finite-dim linear systems."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"Again $\\cK$ is linear and we can define the eigensystem as\n",
"$$\n",
"\\cK\\phi = \\mu\\phi\n",
"$$\n",
"and by definition of $\\cK$, $\\mu\\phi=\\dot\\phi$."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "fragment"
}
},
"source": [
"Similar to the discrete case, the eigensystem also has a lattice property, for $\\phi=\\phi_i\\phi_j$\n",
"$$\n",
"\\cK\\phi = \\dot\\phi = \\dot{\\phi}_i\\phi_j + \\phi_i\\dot{\\phi}_j = \\mu_i\\phi_i\\phi_j + \\mu_j\\phi_i\\phi_j = (\\mu_i + \\mu_j)\\phi_i\\phi_j\n",
"$$\n",
"\n",
"> In fact, $\\lambda=\\exp(\\mu\\Delta t)$, which explains the difference in the lattices of $\\cK_t$ (multiplicative) and $\\cK$ (additive)"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"Taking a further look at $\\dot{\\phi}$, by chain rule,\n",
"$$\n",
"\\dot\\phi(\\vx) = \\nabla\\phi(\\vx)\\cdot\\dot{\\vx} = \\nabla\\phi(\\vx)\\cdot\\vf(\\vx)\n",
"$$\n",
"where $\\nabla$ is gradient w.r.t. $\\vx$.\n",
"\n",
"This leads to an alternative eigenvalue problem in the form of PDE\n",
"$$\n",
"\\mu\\phi(\\vx) = \\nabla\\phi(\\vx)\\cdot\\vf(\\vx)\n",
"$$\n",
"\n",
"> This provides another approach to computing $\\phi$ if $\\vf$ is known"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"## EDMD, Cont'd\n",
"\n",
"If $\\vf$ is unknown, which is most of the case, one could develop a continuous version of EDMD starting from the eigenvalue problem\n",
"$$\n",
"\\mu\\phi(\\vx) = \\nabla\\phi(\\vx)\\cdot\\dot{\\vx}\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"If one introduce the approximation\n",
"$$\n",
"\\phi(\\vx) \\approx \\sum_{i=1}^m g_i(\\vx)u_i = \\vg(\\vx)^T\\vu\n",
"$$\n",
"The eigenvalue problem becomes\n",
"$$\n",
"\\begin{align}\n",
"\\mu\\left(\\sum_{i=1}^m g_i(\\vx)u_i\\right) &= \\nabla\\left( \\sum_{i=1}^m g_i(\\vx)u_i \\right)\\cdot\\dot{\\vx} \\\\\n",
"\\mu\\vg(\\vx)^T\\vu &= \\sum_{i=1}^m \\left( \\nabla g_i(\\vx)\\cdot \\dot{\\vx}\\right) u_i\n",
"\\end{align}\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"From a series of data,\n",
"$$\n",
"\\mu\\begin{bmatrix}\\vg(\\vx_1)^T \\\\ \\vg(\\vx_2)^T \\\\ \\cdots \\\\ \\vg(\\vx_N)^T\\end{bmatrix}\\vu =\n",
"\\begin{bmatrix}\n",
"\\nabla g_1(\\vx)\\cdot \\dot{\\vx}_1 & \\nabla g_2(\\vx)\\cdot \\dot{\\vx}_1 & \\cdots & \\nabla g_m(\\vx)\\cdot \\dot{\\vx}_1 \\\\\n",
"\\nabla g_1(\\vx)\\cdot \\dot{\\vx}_2 & \\nabla g_2(\\vx)\\cdot \\dot{\\vx}_2 & \\cdots & \\nabla g_m(\\vx)\\cdot \\dot{\\vx}_2 \\\\\n",
"\\vdots & \\vdots & \\ddots & \\vdots \\\\\n",
"\\nabla g_1(\\vx)\\cdot \\dot{\\vx}_N & \\nabla g_2(\\vx)\\cdot \\dot{\\vx}_N & \\cdots & \\nabla g_m(\\vx)\\cdot \\dot{\\vx}_N\n",
"\\end{bmatrix}\\vu\n",
"$$\n",
"Or in short\n",
"$$\n",
"\\mu\\vG\\vu = \\vtG\\vu\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"Yet another eigenvalue problem!\n",
"\n",
"But since $\\dot{\\vx}$ is involved, a natural idea is to use SINDy, which also helps selecting the observables (i.e. features). See [1] for more details.\n",
"\n",
"[1] Kaiser, E., et al., Data-driven discovery of Koopman eigenfunctions for control, 2021."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"## A touch on Control\n",
"\n",
"Finally let's look at a nonlinear system with control,\n",
"$$\n",
"\\begin{align}\n",
"\\mbox{System dynamics: }&\\ \\dot{\\vx}=\\vf(\\vx,\\vu) \\\\\n",
"\\mbox{Observation: }&\\ \\vg=\\vg(\\vx,\\vu)\n",
"\\end{align}\n",
"$$\n",
"> This topic is still in active research and there are still some fundamental challenges."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"In general, for a scalar observable, one can find\n",
"$$\n",
"\\dot{g} = \\cK^\\vu g\n",
"$$\n",
"where $\\cK^\\vu$ is a linear non-autonomous Koopman operator, that depends on the control input $\\vu$.\n",
"\n",
"Like before, let's consider the eigenfunctions\n",
"$$\n",
"\\dot\\phi = \\cK^\\vu \\phi = \\mu^\\vu \\phi\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"If we factor in the actual dynamics $\\vf(\\vx,\\vu)$\n",
"$$\n",
"\\mu^\\vu \\phi = \\cK^\\vu \\phi = \\dot{\\phi} = \\nabla_\\vx \\phi \\cdot \\vf(\\vx,\\vu) + \\nabla_\\vu \\phi \\cdot \\dot{\\vu}\n",
"$$\n",
"> Problem: What is $\\dot\\vu$?\n",
"\n",
"This produces [non-causal system](https://en.wikipedia.org/wiki/Causal_system) that is inconvenient to analysis, even in classical control methods.\n",
"\n",
"+ So we set $\\dot\\vu=0$ and assume step changes in control inputs (i.e. zero-hold)"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"One gets\n",
"$$\n",
"\\dot{\\phi} = \\nabla_\\vx \\phi \\cdot \\vf(\\vx,\\vu)\n",
"$$\n",
"\n",
"> Problem: $\\vf(\\vx,\\vu)$ is too generic, can we try something more tangible?\n",
"\n",
"+ So to make life even easier, we assume a control-affine system,\n",
"$$\n",
"\\dot{\\vx} = \\vf_0(\\vx) + \\sum_{i=1}^L \\vf_i(\\vx)u_i\n",
"$$\n",
"where the system is linear in $\\vu$.\n",
"+ The autonomous part $\\vf^0$ has a Koopman form\n",
"$$\n",
"\\mu^\\vf\\phi = \\cK^\\vf \\phi = \\dot{\\phi} = \\nabla_\\vx \\phi \\cdot \\vf_0(\\vx)\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"The eigenproblem simplifies to\n",
"$$\n",
"\\begin{align}\n",
"\\dot\\phi &= \\nabla_\\vx \\phi \\cdot \\vf_0(\\vx) + \\sum_{i=1}^L \\nabla_\\vx \\phi \\cdot \\vf_i(\\vx)u_i \\\\\n",
"&= \\mu^\\vf \\phi + \\sum_{i=1}^L \\nabla_\\vx \\phi \\cdot \\vf_i(\\vx)u_i\n",
"\\end{align}\n",
"$$\n",
"\n",
"> Problem: What about the second part?\n",
"\n",
"+ Let's assume $\\nabla_\\vx \\phi \\cdot \\vf_i(\\vx)$ can be approximated by a specific set of eigenfunctions $\\{\\phi_j\\}_{j=1}^P$\n",
"$$\n",
"\\nabla_\\vx \\phi \\cdot \\vf_i(\\vx) = \\sum_{j=1}^P b_{ij}\\phi_j\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"Then for the $k$th eigenfunction,\n",
"$$\n",
"\\dot{\\phi}_k = \\mu_k^\\vf\\phi_k + \\sum_{i=1}^L \\sum_{j=1}^P b_{ij}^k\\phi_j u_i \\equiv \\mu_k^\\vf\\phi_k + (\\vB^k\\vtf)^T \\vu\n",
"$$\n",
"And if one write the equation for the $P$ eigenfunctions,\n",
"$$\n",
"\\dot\\vtf = \\vtL\\vtf + \\begin{bmatrix} \\vB^1\\vtf & \\vB^2\\vtf & \\cdots & \\vB^P\\vtf \\end{bmatrix}^T\\vu = \\vtL\\vtf + (\\cB\\vtf) \\vu\n",
"$$\n",
"+ $\\vtL$ is a diagonal matrix with eigenvalues $\\{\\mu_j\\}_{j=1}^P$\n",
"+ $\\cB$ is a \"3D matrix\", or a so-called bilinear operator\n",
"\n",
"> Now we get a working **Koopman Bilinear Form**\n",
"\n",
"+ We can apply a modified EDMD technique to find the coefficients and the eigenfunctions, see [Brunton2022]"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "slide"
}
},
"source": [
"## Two examples"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"#### Example 1\n",
"\n",
"\n",
"\n",
"Ref: Li, Y., et al., Learning compositional Koopman operators for model-based control, 2019."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"\n",
"\n",
"Ref: Li, Y., et al., Learning compositional Koopman operators for model-based control, 2019."
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"#### Example 2\n",
"\n",
"\n",
"\n",
"In-house work"
]
},
{
"cell_type": "markdown",
"metadata": {
"slideshow": {
"slide_type": "subslide"
}
},
"source": [
"\n",
"\n",
"In-house work"
]
}
],
"metadata": {
"celltoolbar": "Slideshow",
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.6"
}
},
"nbformat": 4,
"nbformat_minor": 2
}