|
| 1 | +# %% |
| 2 | +%matplotlib notebook |
| 3 | +import numpy as np |
| 4 | +import matplotlib.pyplot as plt |
| 5 | +plt.rcParams.update({ 'text.usetex': False, 'mathtext.fontset': 'cm', |
| 6 | + 'font.size': 12.0, 'axes.labelsize': 'medium', |
| 7 | + 'xtick.labelsize': 'x-small', 'ytick.labelsize': 'x-small', |
| 8 | + 'axes.grid': True, 'axes.formatter.limits': [-3, 6], |
| 9 | + 'grid.alpha': 0.5, 'figure.figsize': [11.0, 4], |
| 10 | + 'figure.constrained_layout.use': True, 'scatter.marker': 'x', |
| 11 | + 'animation.html': 'jshtml'}) |
| 12 | + |
| 13 | +from IPython.display import display, Markdown |
| 14 | + |
| 15 | +import warnings |
| 16 | +warnings.filterwarnings("ignore") |
| 17 | + |
| 18 | +import scipy.signal as signal |
| 19 | +import scipy.linalg as sclin |
| 20 | +import numpy.random as rng |
| 21 | +import numpy.linalg as lin |
| 22 | +import control as cm |
| 23 | +from helperFunctions import * |
| 24 | + |
| 25 | +# %% [markdown] |
| 26 | +# # Time domain |
| 27 | +# ## Systems and differential equations |
| 28 | +# *The basics of all things control* |
| 29 | +# |
| 30 | +# *Systems are all and all is a system.* So you can decribe dynamic processses with physics, then you get equations with time derivatives. So for Linear-Time-Independent (LTI) systems, the differential equation with output $y(t)$ and input $u(t)$ is |
| 31 | +# $$ \frac{d^n}{d^nt}y + a_1\frac{d^{n-1}}{d^{n-1}t}y + \cdots + a_n y = \frac{d^m}{d^mt}u + b_1\frac{d^{m-1}}{d^{m-1}t}u + \cdots + b_m u.$$ |
| 32 | +# Here, $n$ is called the order of the system. |
| 33 | +# |
| 34 | +# We want to be able to solve Initial Value Problems (IVPs), because that's cool and it shows us what the system will do from some initial state. This solution is built in two steps, the homogenous solution and satisfying the initial conditions. |
| 35 | +# |
| 36 | +# The homogenous solution is obtained by equation the left hand side to zero. For this, we'll steal a little ahead from the Laplace transform where |
| 37 | +# $$ \frac{d^n}{d^nt} \overset{\mathcal{F}}{\rightarrow} s^n. $$ |
| 38 | +# Applying this and eliminating $y$ to obtain the homogenous solution yields |
| 39 | +# $$\frac{d^n}{d^nt}y + a_1\frac{d^{n-1}}{d^{n-1}t}y + \cdots + a_n y = 0 $$ |
| 40 | +# $$\overset{\mathcal{F},\; \frac{1}{y}}{\rightarrow} s^n + a_1 s^{n-1} + \cdots + a_n = 0. $$ |
| 41 | +# This is a polynomial of order $n$ called **the characteristic polynomial** and we know then it has $n$ roots as well. |
| 42 | +# |
| 43 | +# These roots, $\lambda_k$, actually form the solution to the homogenous problem, since a polynomial is the product of its roots, |
| 44 | +# $$ \prod_{k=1}^n (s-\lambda_k).$$ |
| 45 | +# The solution to these types of ODEs are sums of exponentials of the form $c_k e^{\lambda_k t}$, where $c_k$ is determined through the initial condition. Therefore, the complete solution is |
| 46 | +# $$ y(t) = \sum_{k= 1}^n c_k e^{\lambda_k t} .$$ |
| 47 | +# |
| 48 | +# If all $\mathfrak{R}(\lambda_k) < 0$ then $y(t)$ goes to 0 and the system is stable. The Routh-Hurwitz criterion gives the stability requirements up to the third degree. Stability is guaranteed for systems of order |
| 49 | +# 1. if $a_1 > 0$. |
| 50 | +# 2. if $a_1, a_2 > 0$. |
| 51 | +# 3. if $a_1, a_2, a_3 >0$ and $ a_1a_2 > a_3$. |
| 52 | +# |
| 53 | +# Real-valued $\lambda_k$ yield an exponential trajectory $e^{\lambda_k t}$. Complex-values come in conjugate pairs so they yield an exponential trajectory multiplied by a cosine. Because $\lambda_k, \lambda_{k+1} = \epsilon \pm j\omega$, the trajectory of this pair becomes $e^{(\epsilon \pm j\omega)t} = e^{\epsilon t} (e^{j\omega t} + e^{-j\omega t}) = e^{\epsilon t} (2\cos(\omega t))$. Play around with the code block here to get a feel for it! |
| 54 | +# |
| 55 | + |
| 56 | +# %% |
| 57 | +## Eigenvalues |
| 58 | +sig1 = -0.8 |
| 59 | +sig2 = 3 + 15j |
| 60 | +sig3 = -0.5 + 5j |
| 61 | + |
| 62 | +###### Plotting ######### |
| 63 | +fig, ax = plt.subplots(1,3) |
| 64 | +for sig, idx in zip([sig1, sig2, sig3], range(3)): |
| 65 | + t = np.linspace(0, 3/abs(np.real(sig)), num=300) # Adapt to convergence speed |
| 66 | + if np.iscomplex(sig): # Plot decomposition |
| 67 | + l1, = ax[idx].plot(t, np.exp(sig.real * t), 'r--', label=f"$e^{r"{"}{sig.real}t{r"}"}$") # Upper envelope |
| 68 | + ax[idx].plot(t, -np.exp(sig.real * t), 'r--') # Lower envelope |
| 69 | + l2, = ax[idx].plot(t, 2*np.cos(sig.imag * t), 'k', alpha=0.2, label=f"$2\cos({sig.imag} t)$") # Oscillation |
| 70 | + ax[idx].legend(handles=[l1, l2]) |
| 71 | + ax[idx].plot(t, np.exp(sig * t), 'k') # Trajectory |
| 72 | + ax[idx].set(title=f"$\lambda = {sig}$", xlabel="$t$") |
| 73 | + |
| 74 | +ax[0].set_ylabel("$y(t)$") |
| 75 | +display(fig) |
| 76 | + |
| 77 | + |
| 78 | +# %% [markdown] |
| 79 | +# ## State space representation |
| 80 | +# *A blessing from above* |
| 81 | +# |
| 82 | +# Cool we know systems now, but it's a bit of an ugly system representation to be honest. Luckily we live in a universe where the state space representation exists and we can recast any system as a vector first order ODE and an output equation. This representation is of the form $$ \dot x = f(x, u)$$ $$ y = h(x, u),$$ |
| 83 | +# where $x$ is the state of the system, $u$ the input, and $y$ the output. The state dynamics are described with $f: \mathbb{R}^n\times\mathbb{R}^p\rightarrow\mathbb{R}^n$ and the output measurements with $h: \mathbb{R}^n\times\mathbb{R}^p\rightarrow\mathbb{R}^q$. |
| 84 | +# |
| 85 | +# What black magic do we perform to get these first order ODEs? Well, suppose you have a second order ODE in $v$, $\ddot v = \dot v + v$, then this is equivalent to the first order ODE $$ \begin{bmatrix}\dot v\\\ddot v\end{bmatrix} = \begin{bmatrix}0&1\\1&1\end{bmatrix}\begin{bmatrix}\dot v\\ v\end{bmatrix}.$$ |
| 86 | +# |
| 87 | +# Last but not least of the amazing aspects of the state space representation: there are many nice numerical integrators to simulate them. Think forward Euler or Runge-Kutta. |
| 88 | + |
| 89 | +# %% [markdown] |
| 90 | +# <div style="text-align:center;background-color:tomato;">End of lecture 2 (Lecture 1 was course logistics)</div> |
| 91 | + |
| 92 | +# %% [markdown] |
| 93 | +# ### Block diagrams |
| 94 | +# *Pixel perfect ways to visualise systems* |
| 95 | +# |
| 96 | +#  |
| 97 | +# |
| 98 | +# It's a nice sanity check that any fundamental block scheme of an $n$-th order system has $n$ integrators. |
| 99 | +# |
| 100 | +# |
| 101 | +# ## Equilibrium points |
| 102 | +# *Forever and unchanging* |
| 103 | +# |
| 104 | +# Often, systems will have equilibria, meaning there are states where system will remain over time. These state-input point pairs are denoted as $(x_e, u_e)$ or $(\bar x, \bar u)$. I'll use the latter since this is more traditional in control and we have to appease our control elders (JW). So equilibria remain unchanging, so the derivative of the state is zero. In mathematicians' language |
| 105 | +# $$ \dot{\bar x} = f(\bar x, \bar u) = 0. $$ |
| 106 | +# Finding the equilibria states as a function of the equilibria inputs is as simple as solving this equation. |
| 107 | +# |
| 108 | +# |
| 109 | +# |
| 110 | +# |
| 111 | + |
| 112 | +# %% [markdown] |
| 113 | +# <div style="text-align:center;background-color:tomato;">End of lecture 3</div> |
0 commit comments