Chapter 4 The step response
-
Remark 4.1. We recall that by the variation of parameters formula we have that (3.1) has the unique solution
\[ x(t)=\e ^{At}x^0+\int _0^t \e ^{A(t-\theta )}Bu(\theta )\,d\theta ,\qquad y(t)=C\e ^{At}x^0+\int _0^t C\e ^{A(t-\theta )}Bu(\theta )\,d\theta +Du(t), \]
where \(z\mapsto \e ^{Az}\) is the matrix exponential function defined by
\[ \e ^{Az}=\sum _{k=0}^\infty \frac {1}{k!}A^kz^k. \]
-
Definition 4.2. The step response of (3.1) is the \(\mR ^{p\times m}\) valued function
\[ H(t)=\begin {dcases} 0&t<0\\ \int _0^t C\e ^{A\theta }B\,d\theta +D&t\geq 0. \end {dcases} \]
-
Remark 4.3. Assume that \(m=p=1\) (the single-input-single-output (SISO) case). By the variation of parameters formula the output for the input (and initial condition zero)
\(\seteqnumber{0}{4.}{0}\)\begin{equation} \label {eq:step} u(t)=\begin{dcases} 0&t<0\\ 1&t\geq 0, \end {dcases} \end{equation}
(this \(u\) is called the unit step or Heaviside function) equals the step response \(H\). This is the origin of the name step response. By arguing in terms of components, this can be extended to the case of general \(m\) and \(p\) (the multi-input-multi-output (MIMO) case).
-
Remark 4.4. We will see later that the step response uniquely determines the output for all inputs (when \(x^0=0\)), so it completely characterizes the input-output behavior of the system (3.1).
4.1 Examples
4.1.1 First order systems
Consider the first order system
\[ T\dot {y}(t)+y(t)=Ku(t), \]
where \(T,K>0\). The parameter \(T\) is called the time constant and \(K\) the steady state gain of the system. By defining \(x:=y\), this can be written as an input-state-output system with \(n=m=p=1\), \(A=-\frac {1}{T}\), \(B=\frac {K}{T}\), \(C=1\) and \(D=0\). The step response of this system is (for \(t>0\))
\[ H(t)=K\left (1-\e ^{-t/T}\right ). \]
If we re-scale time by defining \(\theta :=t/T\) then we obtain (using \(\frac {d}{dt}=\frac {d\theta }{dt}\frac {d}{d\theta }\))
\[ T\frac {1}{T}\frac {dy}{d\theta }+y=Ku, \]
i.e.
\[ \frac {dy}{d\theta }+y=Ku. \]
If we then re-scale \(y\) by defining \(\tilde {y}:=y/K\) then we obtain
\[ \frac {d\tilde {y}}{d\theta }+\tilde {y}=u. \]
The step response of this system is
\[ \tilde {H}(\theta )=1-\e ^{-\theta }, \]
so that up to scaling the axes, the step response of a physical first order system is \(1-\e ^{-t}\). The step response for a first order system is given in Figure 4.1. In Figures 4.2a, 4.2b we illustrate how the step response depends on the steady state gain \(K\) and the time constant \(T\).
4.1.2 Case study: diagnosis of lung diseases
One method of diagnosing lung diseases (and distinguishing between two common types of lung diseases) is spirometry (see Figure 4.3a). There a patient blows as hard as possible for 6 to 8 second and the volume exhaled is measured. That the patient blows as hard as possible is not actually so crucial, more important is that the blowing is at a constant level (and the easiest way to make sure of this is to ask the patient to blow as hard as possible). The blowing acts as an input and the volume exhaled is the output. Since the blowing is at a constant level, the input is a step (or at least: the first 6 to 8 seconds of a step) and therefore the output is the step response (or at least: the first 6 to 8 seconds of the step response). The underlying input-state-output system is the lungs. Therefore what spirometry measures is the step response of the lungs of the patient. Note that by Remark 4.4 the step response uniquely characterizes the input-output behavior, so even though spirometry seems restricted to one very specific input (blowing as hard as possible for 6 to 8 seconds), it actually captures the entire input-output behavior of the lungs.
From Figure 4.3b we see that the lungs are essentially a first order system (comparing the shapes of the graphs to Figure 4.1). Comparing Figure 4.3b to Figure 4.2a we see that for restrictive lung diseases the steady state gain is lower than normal (whereas the time constant is the same as normal). In Figure 4.3b we also see that for obstructive lung diseases (such as asthma) the steady state gain is lower than normal; however the literature is not unanimous about this. Comparing Figure 4.3b to Figure 4.2b we see that for obstructive lung diseases (such as asthma) the time constant is larger than normal. This is consistent in the literature. So the combination of gain and time constant distinguishes normal lungs and two types of lung diseases.
In Figure 4.3b FVC stands for forced vital capacity and is what in general is called the steady state gain. FEV1 stands for forced expiratory volume in one second. The ratio FEV1/FVC is equivalent to the time constant (in the sense that one can be uniquely computed from the other). What is more commonly measured is the ratio FEV1/FEV6 which is also equivalent to the time constant (the important thing is that we look at a ratio so that the steady state gain cancels).
4.1.3 Second order systems (low pass)
Consider the second order system
\(\seteqnumber{0}{4.}{1}\)\begin{equation} \label {eq:secondlow} T^2\ddot {y}(t)+2\zeta T\dot {y}(t)+y(t)=Ku(t), \end{equation}
where \(T,K>0\) and \(\zeta \geq 0\). The parameter \(T\) is called the time constant and \(K\) the steady state gain. The parameter \(\zeta \) is the damping ratio. The equation can alternatively by written as
\[ \ddot {y}(t)+2\zeta \omega _0\dot {y}(t)+\omega _0^2y(t)=K\omega _0^2u(t), \]
where \(\omega _0=\frac {1}{T}\) is called the natural frequency. By defining
\[ n=2,\quad m=p=1,\qquad x=\bbm {y\\\dot {y}},\quad A=\bbm {0&1\\-\omega _0^2&-2\zeta \omega _0},\quad B=\bbm {0\\K\omega _0^2},\quad C=\bbm {1&0},\quad D=0. \]
this can be written as an input-state-output system. If we re-scale \(\theta :=t/T\) and \(\tilde {y}:=y/K\), then we obtain
\[ \frac {d^2\tilde {y}}{d\theta ^2}+2\zeta \frac {d\tilde {y}}{d\theta }+\tilde {y}=u. \]
In contrast to the first order case, the parameter \(\zeta \) still appears in the re-scaled differential equation and the form of the step response crucially depend on \(\zeta \).
The roots of \(s^2+2\zeta s+1\) as a function of \(\zeta \) are given in Figure 4.4. If \(\zeta \in [0,1)\phantom {]}\), then \(s^2+2\zeta s+1\) has two complex conjugate (non-real) roots. This is called the underdamped case. If \(\zeta =1\), then \(s^2+2\zeta s+1\) has a double (real) root. This is called the critically damped case. If \(\zeta >1\), then \(s^2+2\zeta s+1\) has two distinct real roots. This is called the overdamped case.
The step responses for relevant values of \(\zeta \) are given in Figure 4.5. When \(\zeta \geq 1\), the step response never crosses the horizontal line at height \(K\). Instead when \(\zeta \in [0,1)\phantom {]}\), the step response does cross this line, i.e. the step response is oscillatory.
-
Example 4.5. We calculate the step response of
\[ \ddot {y}+\dot {y}+y=u. \]
We only need to calculate for \(t>0\) (as the step response is zero for \(t<0\)), so we solve
\[ \ddot {y}+\dot {y}+y=1,\qquad y(0)=0,\quad \dot {y}(0)=0. \]
A particular solution is \(y=1\). The homogeneous equation is
\[ \ddot {y}+\dot {y}+y=0. \]
The characteristic polynomial of this is \(s^2+s+1\), which has roots \(\frac {-1\pm i\sqrt {3}}{2}\), so that the general solution of the homogeneous equation is
\[ \e ^{-t/2}\left (a \cos \left (\frac {\sqrt {3}}{2}t\right )+b \sin \left (\frac {\sqrt {3}}{2}t\right )\right ). \]
adding the particular solution gives the general solution of our original equation:
\[ y(t)=1+\e ^{-t/2}\left (a \cos \left (\frac {\sqrt {3}}{2}t\right )+b \sin \left (\frac {\sqrt {3}}{2}t\right )\right ). \]
We now use the (zero) initial conditions to determined the constants \(a\) and \(b\). We have
\[ 0=y(0)=1+a,\quad 0=\dot {y}(0)=\frac {-1}{2}a+\frac {\sqrt {3}}{2}b. \]
This gives \(a=-1\) and \(b=\frac {-1}{\sqrt {3}}\). Therefore the step response is
\[ H(t)=1+\e ^{-t/2}\left (-\cos \left (\frac {\sqrt {3}}{2}t\right )-\frac {1}{\sqrt {3}} \sin \left (\frac {\sqrt {3}}{2}t\right )\right ). \]
The formulas for the step response in general are as follows:
\(\seteqnumber{0}{4.}{2}\)\begin{gather*} H(t)=K\left [1-\e ^{-\zeta t/T}\left (\cosh (\sqrt {\zeta ^2-1}~(t/T))+\frac {\zeta }{\sqrt {\zeta ^2-1}}\sinh (\sqrt {\zeta ^2-1}~(t/T))\right )\right ] \text {when~} \zeta >1, \\ H(t)=K\left [1-\e ^{-t/T}(1+t/T)\right ] \text {when~} \zeta =1, \\ H(t)=K\left [1-\e ^{-\zeta t/T}\left (\cos (\sqrt {1-\zeta ^2}(t/T))+\frac {\zeta }{\sqrt {1-\zeta ^2}}\sin (\sqrt {1-\zeta ^2}(t/T))\right )\right ] \text {when~} \zeta \in [0,1)\phantom {]}. \end{gather*}
-
Example 4.6. Consider a mass-spring-damper system with a force input and position output:
\[ m\ddot {y}+d\dot {y}+ky=u. \]
We can write this in the standard form (4.2) where
\[ K=\frac {1}{k},\qquad T=\sqrt {\frac {m}{k}},\qquad \zeta =\frac {d}{2\sqrt {km}}. \]
Consider the situation where the mass is suspended from the (fixed) ceiling to which it is attached by the spring and damper. Assume that at time zero the system is supported from below to be at its equilibrium position. At time zero this support is removed and the constant gravitational force \(u=mg\) now acts on it. The output then is the step response multiplied by \(mg\). The new equilibrium position then is the steady state gain times this constant force: \(\frac {mg}{k}\) (we measure \(y\) as being positive towards the ground). If the damping constant is large, so that \(\zeta >1\) (the overdamped case), then convergence to this new equilibrium value is monotonic. If the damping constant is small (so that \(\zeta \in (0,1)\)), then the mass overshoots the new equilibrium on its approach.
4.1.4 Second order systems (band pass)
Consider the second order system
\(\seteqnumber{0}{4.}{2}\)\begin{equation} \label {eq:secondband} T^2\ddot {q}(t)+2\zeta T\dot {q}(t)+q(t)=Ku(t),\qquad y(t)=\dot {q}(t), \end{equation}
where again \(T,K>0\) and \(\zeta \geq 0\). By defining
\[ n=2,\quad m=p=1,\qquad x=\bbm {q\\\dot {q}},\quad A=\bbm {0&1\\-\omega _0^2&-2\zeta \omega _0},\quad B=\bbm {0\\K\omega _0^2},\quad C=\bbm {0&1}, \quad D=0, \]
this can be written as an input-state-output system (note that the only difference with Section 4.1.3 is in \(C\)).
The step response of (4.3) is the derivative of the step response of 4.2. The step responses for relevant values of \(\zeta \) are given in Figure 4.6.
The step response of (4.3) in the overdamped case for varying \(K\) is given in Figure 4.7a and for varying \(T\) is given in Figure 4.7b.
4.1.5 Case study: muscle energy
To perform work, the muscles need energy. This can be provided through the aerobic energy system (which needs oxygen) and through two anaerobic energy systems (which do not need oxygen). To investigate these energy systems, experimenters put someone on a stationary bike and ask them to put in a maximal effort for about 100 seconds. The cyclist is attached to various devices for measurement purposes. That the effort is maximal is not crucial, important is that the effort is constant (and an easy way to achieve this is to ask for maximal effort). The input then is a step and therefore the measurements give the step response of the three energy systems. These are given in Figure 4.8a.
The aerobic energy system (labelled “mitochindrial respiration” in Figure 4.8a) is a first order system (as we essentially already saw in a different context in Section 4.1.2). The two anaerobic energy systems (labelled “phosphagen” and “glycolytic” in Figure 4.8a) have step responses of overdamped second order systems of the form (4.3). Comparing to Figure 4.7b, we see that the phosphagen (or alactic) energy system has a smaller time constant than the glycolytic (or lactic) energy system.
That the glycolytic (or lactic) energy system is an overdamped second order systems of the form (4.3) can be explained as follows. Denote the power produced by the lactic system by \(P\). The differential equation for \(P\) then is
\[ T_1\dot {P}+P+L=K_1u. \]
This is a standard first order differential equation with time constant \(T_1\) and gain \(K_1\), except for the presence of term \(L\). This term represents the lactic acid in the muscles which inhibits power production by the lactic energy system. Lactic acid accumulates due to earlier power produced:
\[ T_2\dot {L}=P. \]
Differentiating the equation for \(L\) and substituting the equation for \(P\) gives
\[ T_1T_2\ddot {L}+T_2\dot {L}+L=K_1u,\qquad P=T_2\dot {L}. \]
Defining \(q:=T_2L\), this is a second order system of the form (4.3) with
\[ K=T_2K_1,\qquad T=\sqrt {T_1T_2},\qquad \zeta =\sqrt {\frac {T_2}{4T_1}}. \]
This is an overdamped system if and only if \(T_2>4T_1\).
Figure 4.8b gives a comparison between sprinters and endurance athletes. We see that for the anaerobic systems the time constants seem similar, but the gain for the sprinters is higher. For the aerobic system (which is first order) the time constant for sprinters is smaller, but the steady state gain is lower. For the overall energy (labelled “energy demand” in Figure 4.8b) this results in sprinters being able to generate more power for the first 30 seconds, but less power after that.
4.2 Problems
Consider as in Section 1.4 (with \(m=1\), \(d=0\) and \(k=1\))
\[ \dot {q}=v-v_e,\qquad \dot {v}+q=F_e. \]
In the context of this chapter we have the input \(u:=\sbm {F_e\\v_e}\) and the output \(y:=\sbm {q\\v}\) (this is therefore a multi-input-mult-output system). Determine the step response of the above system.
-
Solution. We only have to consider \(t>0\) (since the step response is zero for \(t<0\)). We determine the step response matrix column-by-column. We first consider \(F_e=1\) and \(v_e=0\) (and zero initial condition). This gives
\[ \dot {q}=v,\qquad \dot {v}+q=1,\qquad q(0)=0,\quad v(0)=0. \]
We can rewrite this as the second order equation (we could alternatively solve the system of equations using matrix methods):
\[ \ddot {q}+q=1,\qquad q(0)=0,\quad \dot {q}(0)=0. \]
a particular solution is \(q=1\). The homogeneous equation is \(\ddot {q}+q=0\) and has general solution \(a\cos (t)+b\sin (t)\). Therefore the general solution of our original equation is
\[ 1+a\cos (t)+b\sin (t). \]
From the zero initial conditions we obtain \(1+a=0\), \(b=0\). Therefore we obtain
\[ q(t)=1-\cos (t). \]
We then calculate \(v=\dot {q}=\sin (t)\). Therefore the first column of the step response matrix is
\[ \bbm {1-\cos (t)\\\sin (t)}. \]
To obtain the second column of the step response matrix, we consider \(F_e=0\) and \(v_e=1\) (and zero initial conditions). This gives
\[ \dot {q}=v-1,\qquad \dot {v}=-q,\qquad q(0)=0,\quad v(0)=0. \]
a particular solution is \(v=1\), \(q=0\). The homogeneous equation is
\[ \dot {q}=v,\qquad \dot {v}=-q, \]
which we can write as the second order system
\[ \ddot {q}=-q, \]
the general solution of which is \(q(t)=a\cos (t)+b\sin (t)\). This gives \(v(t)=\dot {q}(t)=-a\sin (t)+b\cos (t)\). Therefore the general solution of the original equation is
\[ q(t)=a\cos (t)+b\sin (t),\qquad v(t)=1-a\sin (t)+b\cos (t). \]
We now determine \(a\) and \(b\) from the zero initial conditions. This gives \(a=0\), \(1+b=0\). Therefore
\[ q(t)=-\sin (t),\qquad v(t)=1-\cos (t). \]
It follows that the second column of the step response is
\[ \bbm {-\sin (t)\\1-\cos (t)}. \]
Hence the step response is
\[ H(t)=\bbm {1-\cos (t)&-\sin (t)\\\sin (t)&1-\cos (t)}. \]
We can alternatively use the formula from Definition 4.2. We write our system of equations in matrix form, with \(x:=\sbm {q\\v}\), as
\[ \dot {x}=\bbm {0&1\\-1&0}x+\bbm {0&-1\\1&0}u,\qquad y=x. \]
Therefore
\[ A=\bbm {0&1\\-1&0},\qquad B=\bbm {0&-1\\1&0},\qquad C=\bbm {1&0\\0&1},\qquad D=\bbm {0&0\\0&0}. \]
The matrix exponential function of \(A\) is (the second year module MA20220 Ordinary Differential Equations and Control gave several different ways of calculating this):
\[ \e ^{At}=\bbm {\cos (t)&\sin (t)\\-\sin (t)&\cos (t)}. \]
We therefore have
\(\seteqnumber{0}{4.}{3}\)\begin{align*} H(t)&=\int _0^t C\e ^{A\theta }B\,d\theta +D =\int _0^t\bbm {\cos (\theta )&\sin (\theta )\\-\sin (\theta )&\cos (\theta )}\bbm {0&-1\\1&0}\,d\theta =\int _0^t \bbm {\sin (\theta )&-\cos (\theta )\\\cos (\theta )&\sin (\theta )}\,d\theta \\&=\bbm {1-\cos (t)&-\sin (t)\\\sin (t)&1-\cos (t)}, \end{align*} which is the same as we computed above in a column-by-column fashion. □