Time Delay in Dynamic Systems
Time delay is a shift in the effect of an input on an output dynamic response. A first-order linear system with time delay is:
$$\tau_p \frac{dy(t)}{dt} = -y(t) + K_p u\left(t-\theta_p\right)$$
has variables y(t) and u(t) and three unknown parameters with
$$K_p \quad \mathrm{= Process \; gain}$$
$$\tau_p \quad \mathrm{= Process \; time \; constant}$$
$$\theta_p \quad \mathrm{= Process \; dead \; time}$$
The time delay `\theta_p` is expressed as a time shift in the input variable u(t).
$$u\left(t-\theta_p\right)$$
$$y(t \lt \theta_p) = y(0)$$
$$y(t \ge \theta_p) = \left( e^{-\left(t - \theta_p \right) / \tau_p}\right) y(0) + \left( 1 - e^{-\left(t - \theta_p \right) / \tau_p} \right) K_p \Delta u $$
For a step change `\Delta u`, the analytical solution for a first-order linear system without time delay (x(t) = y(t) with `theta_p=0`)
$$\tau_p \frac{dx(t)}{dt} = -y(t) + K_p u\left(t\right)$$
is
$$x(t) = K_p \left(1-\exp\left(\frac{-t}{\tau_p}\right)\right) \Delta u$$
With dead-time, the solution becomes:
$$y(t) = x \left( t-\theta_p \right) S \left( t-\theta_p \right)$$
$$\quad = K_p \left( 1- \exp\left(-\frac{ t-\theta_p } { \tau_p } \right) \right) \, \Delta u \, S \left( t-\theta_p \right) $$
where
$$S\left(t-\theta_p\right)$$
is a step function that changes from zero to one at `t=theta_p`.
Analytic Solution Derivation with Laplace Transforms
Start with the linear differential equation with time delay:
$$\tau_p \frac{dy(t)}{dt} = -y(t) + K_p u\left(t-\theta_p\right)$$
Perform a Laplace transform from the tables on each part of the equation:
$$\mathcal{L}\left(\tau_p \frac{dy(t)}{dt}\right) = \tau_p \left(s Y(s) - y(0)\right)$$
$$\mathcal{L}\left(-y(t)\right) = -Y(s)$$
$$\mathcal{L}\left(K_p u\left(t-\theta_p\right)\right) = K_p U(s) e^{-\theta_p s}$$
If the input `U(s)` is a step function of size `\Delta u` then:
$$U(s) = \frac{\Delta u}{s}$$
Combining all of the individual Laplace transforms, the equation in Laplace domain with zero initial condition y(0)=0 is then:
$$\tau_p \, s \, Y(s) = -Y(s) + K_p \frac{\Delta u}{s} e^{-\theta_p s}$$
At this point, the equation is algebraic and can be solved for Y(s) by combining terms on each side:
$$\tau_p \, s \, Y(s) + Y(s) = K_p \frac{\Delta u}{s} e^{-\theta_p s}$$
and factoring out the Y(s) term:
$$Y(s) \left(\tau_p \, s + 1\right) = K_p \frac{\Delta u}{s} e^{-\theta_p s}$$
A final steps are to isolate Y(s) on the left side of the equation and perform an inverse Laplace transform to return to the time domain:
$$Y(s) = K_p \frac{\Delta u}{s\,\left(\tau_p \, s + 1\right)} e^{-\theta_p s}$$
$$\mathcal{L}^{-1}\left(Y(s)\right) = y(t) = K_p \left( 1-\exp \left( -\frac{ \left( t-\theta_p \right) } { \tau_p } \right) \right) \, \Delta u \, S \left( t-\theta_p \right)$$
import matplotlib.pyplot as plt
from scipy.integrate import odeint
# specify number of steps
ns = 100
# define time points
t = np.linspace(0,ns/10.0,ns+1)
class model(object):
# process model
Kp = 2.0
taup = 200.0
thetap = 0.0
def process(y,t,u,Kp,taup):
# Kp = process gain
# taup = process time constant
dydt = -y/taup + Kp/taup * u
return dydt
def calc_response(t,m):
# t = time points
# m = process model
Kp = m.Kp
taup = m.taup
thetap = m.thetap
# specify number of steps
ns = len(t)-1
delta_t = t[1]-t[0]
# storage for recording values
op = np.zeros(ns+1) # controller output
pv = np.zeros(ns+1) # process variable
# step input
op[10:]=2.0
# Simulate time delay
ndelay = int(np.ceil(thetap / delta_t))
# loop through time steps
for i in range(0,ns):
# implement time delay
iop = max(0,i-ndelay)
y = odeint(process,pv[i],[0,delta_t],args=(op[iop],Kp,taup))
pv[i+1] = y[-1]
return (pv,op)
# calculate step response
model.Kp = 2.5
model.taup = 2.0
model.thetap = 3.0
(pv,op) = calc_response(t,model)
pv2 = np.zeros(len(t))
for i in range(len(t)):
pv2[i] = model.Kp * (1.0 - np.exp(-(t[i]-model.thetap-1.0)/model.taup))*2.0
pv3 = np.zeros(len(t))
for i in range(len(t)):
pv3[i] = model.Kp * (1.0 - np.exp(-(t[i]-1.0)/model.taup))*2.0
plt.figure(1)
plt.subplot(2,1,1)
frame1 = plt.gca()
frame1.axes.get_xaxis().set_visible(False)
plt.plot(t,pv3,'r--',linewidth=1,label=r'$y(t)$')
plt.plot(t,pv2,'k.-',linewidth=2,label=r'$y(t-\theta_p)$')
plt.plot([0,4,4.0001,10],[0,0,1,1],'g:',linewidth=3,label=r'$S(t- \theta _p)$')
plt.plot(t,pv,'b-',linewidth=4,label=r'$x(t)$')
plt.legend(loc='best')
plt.ylabel('Process Output')
plt.ylim([-4,5])
plt.subplot(2,1,2)
frame1 = plt.gca()
frame1.axes.get_xaxis().set_visible(False)
plt.plot(t,op,'r-',linewidth=3,label=r'$u(t)$')
plt.plot(t+3.0,op,'k--',linewidth=3,label=r'$u(t-\theta_p)$')
plt.ylim([-0.1,2.1])
plt.xlim([0,10])
plt.legend(loc='best')
plt.ylabel('Process Input')
plt.xlabel('Time')
plt.savefig('output.png')
plt.show()