7. Dynamic Systems - Part 1

import numpy as np
import pandas as pd
import scipy as sp
import scipy.signal as sg
from scipy.integrate import odeint
import matplotlib.pyplot as plt
import plotly.express as px
import plotly.graph_objects as go
from plotly.subplots import make_subplots
from ipywidgets import widgets, interact
from IPython import display
---------------------------------------------------------------------------
ModuleNotFoundError                       Traceback (most recent call last)
Input In [1], in <cell line: 2>()
      1 import numpy as np
----> 2 import pandas as pd
      3 import scipy as sp
      4 import scipy.signal as sg

ModuleNotFoundError: No module named 'pandas'

Let us consider three scenarios in which we wish to measure the deflection of a beam -

  • Static response - deflection under constant load

  • Quasi-static response - the deflection as a response to a slowly increasing load, there is enough time for the beam to go from one equilibrium state to the other.

  • Dynamic response - a vibrating beam

In general, when the transient characteristics become important the measurment will be become more complex

You will find, that for different systems designs will respod in different manner to the various types of input signals., thus dictating their suitability for the measurment you are planning. The mathematical modeling of a measurment system is then crucial to choose the proper one.

Reconstructing the input signal of a dynamic signal, we will look for the signal frequency, amplitude and in general the waveform.

It is important to always keep in mind what is it that we are trying to measure and how it interacts with the sensor we chose.

For example, measuring your body heat (~constant) is done through the introduction of a step function to the measuring system….

N\(^{th}\) order systems

An \(N^{th}\) order system with an output signal \(y(t)\), when subjected to a general input signal represented by the forcing signal \(F(t)\) is given by:

\[\begin{split} a_n \frac{d^n y}{dt^n}+a_{n-1} \frac{d^{n-1} y}{dt^{n-1}}+\cdot \cdot \cdot +a_{1} \frac{d y}{dt} + a_0 y= F(t) \\ F(t) = b_m \frac{d^m x}{dt^m}+b_{m-1} \frac{d^{m-1} x}{dt^{m-1}}+\cdot \cdot \cdot +b_{1} \frac{d x}{dt} + b_0 x \quad m \leq n \end{split}\]

and the coefficients \(a_i\) and \(b_j\) are a representation of physical parameters of the measuring system

7.1. N=0

\[ a_0 y(t) = F(t) \]

the constant \(1/a_0\) is called the static sensitivity. In a zeroth order system the system will track the forcing signal instantly.

7.2. N=1

\[\begin{split} a_1 \frac{dy}{dt} +a_0 y(t) = F(t) \\ \tau \frac{dy}{dt} +y(t) =K F(t) \end{split}\]

Solving for a step input

\[\begin{split} F(t) = 0 \quad at \quad t=0 \\ F(t) = A \quad at \quad t \gt 0 \end{split}\]

With the initial condition:

\[ y(t=0)=y_0 \]

we obtain

\[ y(t) = \frac{A}{a_0} + \left( y_0 - \frac{A}{a_0} \right) \exp ^{-t/\tau} \]

7.2.1. dynamic response of a temperature sensor

Tsensor

\(m\) stands for the mass of the liquid in the thermostat, \(c_v\) is the heat capacity (specific) of the liquid, h is related to the convection and \(A\) is the surface area.

We are trying to understand how the liquid’s temperature \(T\) evolves when the thermostat is submerged in an enviroment with temperature \(T_{\infty}\)

From energy conservation:

\[ \frac{dE}{dt} = \dot{Q} = mc_v \frac{dT}{dt} \]

Also, we know that the rate of heat exchange is related to the temprature difference:

\[ \dot{Q} = hA\Delta T \]

and finally

\[ mc_v = hA(T_{\infty}-T) \]

the later translates into a 1\(^st\) order system

\[ \tau \dot{T} +T = T_{\infty} \]

with

\[ \tau = \frac{mc_v}{hA} \]

We wish to understand the response of this system to a step input such that

\[\begin{split} T_{\infty} =\begin{cases} T_1 \ \ , \ \ t<0 \\ T_2 \ \ , \ \ t\geq 0 \end{cases} \end{split}\]

The solution of the above ode reads:

\[ T(t) = T_2 + (T_1 -T_2) \exp \left (\frac{-t}{\tau} \right ) \]

Dynamic error :

\[ e(t) = T(t) - T_2 \]

and when normalized to the initial temperature difference:

\[ \Gamma(t) = \exp \left (\frac{-t}{\tau} \right ) \]

Lets see the effect of changing the time constant of the system again, we will use the package Scipy.signal and pecifically we will use Scipy.signal.lty which is a class for continuous-time linear time invariant systems

%matplotlib inline
def first_order_w_step(tau,totalTime):
    Tvec = np.arange(0.0, totalTime, 0.05)
    FoS = sg.lti(1,[1,1.0/tau])
    step_response = FoS.step(T=Tvec)[1]
    fig, ax = plt.subplots(1, 2, figsize=(18, 8),sharex=True, sharey=True)
    ax[0].plot(Tvec, step_response/step_response[-1],label='system response',color='blue')
    ax[0].axhline(1.0, color='black', label='steady state', linestyle='dashed')
    ax[0].axhline(0.63, color='red',linestyle='dashed' )
    ax[0].axvline(1,color='red')
    ax[0].set_xlabel('Normalized time')
#    ax[0].set_xticks(np.arange(0, totalTime/tau, step=1))
    ax[0].set_title('$1st$ order system response to a step function', fontsize=14)
    ax[0].legend()
    ax[1].plot(Tvec, 1-step_response/step_response[-1],label='Error',color='blue')
    ax[1].axhline(0.0, color='black', label='steady state', linestyle='dashed')
    ax[1].axhline(0.368, color='red',linestyle='dashed' )
    ax[1].axvline(1,color='red')
    ax[1].set_xlabel('Normalized time')
#    ax[1].set_xticks(np.arange(0, totalTime/tau, step=1))
    ax[1].set_title('Error in $1st$ order system response to a step function', fontsize=14)
    ax[1].legend()
    ax[0].set_xlim([0,totalTime/2])
    ax[1].set_xlim([0,totalTime/2])

#totalTime = 2500 #seconds
#tau = 180 #in sec
interact(first_order_w_step, tau = (1,501,10), totalTime = (100,5000,100))
<function __main__.first_order_w_step(tau, totalTime)>

For a sin input:

\[\begin{split} T_{\infty}(t) = T_0 + \Delta T \sin \omega t \\ T(t) = T_0 + B(\omega) \Delta T \sin \left ( \omega t + \phi (\omega) \right ) \\ B(\omega) = (1+(\omega \tau )^2 )^{-1/2} \\ \phi (\omega) = -arctan(\omega \tau) \end{split}\]

In the above set of equations, \(B(\omega)\) and \(\phi(\omega)\) are responsible for controling the amplitude ratio and phase shift between the input and output signals (recall the concept of transfer function).

%matplotlib inline
plt.rcParams['figure.figsize'] = [20, 10]
def first_order_w_sin(tau,Omega):
    Tvec = np.arange(0.0, 2*np.pi/Omega, 0.05)
    tauOmega = np.arange(0.001,1000,0.01)
    FoS = sg.TransferFunction(1,[tau,1.0])
    Ft = np.sin(Omega*Tvec)
    phi=-np.arctan(tauOmega)
    B= np.power(1+np.power(tauOmega,2),-0.5)
    Tout, sinRes, xout = sg.lsim(FoS, Ft, Tvec)
#    circle1 = plt.Circle((tau*Omega,np.power(1+np.power(tau*Omega,2),-0.5)),0.02,color='red')
    fig, ax = plt.subplots(1, 2, figsize=(18, 8),sharex=False, sharey=False)
    ax[0].plot(Tvec*Omega/(2*np.pi), sinRes,label='system response',color='blue')
    ax[0].plot(Tvec*Omega/(2*np.pi), Ft,label='system input',color='green')
    ax[0].set_xlabel('Normalized time')
    ax[0].set_title(r'response to s sin input $\omega \tau$ = %f'%(Omega*tau), fontsize=14)
    ax[0].legend()
    ax[1].plot(tauOmega, B,'-b' )
    ax[1].plot(tau*Omega,np.power(1+np.power(tau*Omega,2),-0.5), 'or', markersize=12)
#    ax[1].add_artist(circle1)
    ax[1].set_ylim(0,1.2)
    ax[1].grid('on')
    ax[1].set_xlabel(r'$\omega \tau$', fontsize=12)
    ax[1].set_xscale('log')

# Simulate tau * dy/dt = -y + F(omega)

interact(first_order_w_sin, tau = (0.1,5,0.1), Omega = (0.1,5,0.1))
<function __main__.first_order_w_sin(tau, Omega)>