Control Analysis and Design through Frequency Domain Methods

Transfer Function Shaping through Control: Closed-Loop vs. Open-Loop

In the relevant section, we discussed some control theoretic configurations mapping an input to an output and how the map can be shaped by control design. Two common architectures are depicted in Figure 1.2 (as a general output feedback; here the control depends on the output of the system) and in Figure 1.3 (as an error output feedback control system; here the control depends on the error between the external input and the system output). More general configurations are also possible, as discussed in the relevant section. For consistency, we will focus on a particular architecture noting that the analysis to follow can be generalized to any of these models.

Consider the (error) feedback loop given in Figure 8.4. Here $P(s)$ denotes the transfer function of the system to be controlled. By writing $Y(s) = P(s)C(s)(R(s) - Y(s))$, it follows that

$\frac{Y(s)}{R(s)} = \frac{P(s)C(s)}{1 + P(s)C(s)}$

is the closed-loop transfer function (under negative unity feedback). Compare this with the setup where there is no feedback: in this case, the transfer function would have been $P(s)C(s)$. This latter expression is often called the (open) loop-transfer function.

The goal is to shape the closed-loop transfer function via the control characterized with $C(s)$ in the frequency domain.

Some motivation via a common class of controllers: PID controllers

By Laplace transforms, we know that differentiation in time domain entails a multiplication with $s$, and integration involves multiplication with $\frac{1}{s}$. In the context of the setup of Figure 8.4, let $e(t) = y(t) - r(t)$. A popular and practical type of control structure involves:

$u(t) = k_i \int_0^t e(t)\, dt + k_d \frac{de}{dt} + k_p e(t)$

which writes as, in Laplace domain,

$U(s) = \left(\frac{k_i}{s} + k_d s + k_p\right) E(s).$

Thus, the control uses both the error itself (proportional), its integration, and its derivative; leading to the term PID control.

---

Bode-Plot Analysis

Bode plots were studied earlier in class. With Bode plots, we observed that we can identify the transfer function of a system, when a system is already stable. However, the Bode plot does not provide insights on how to design a control system or how to adjust a controller so that stability is attained.

---

The Root Locus Method

The set $\{s : 1 + C(s)P(s) = 0\}$ consists of the poles of the transfer function. If one associates with the controller a gain parameter $K$, the root locus method traces the set of all poles as $K$ ranges from $0$ to $\infty$ (and often from $0$ to $-\infty$ as well), so that $\{s : 1 + KC(s)P(s) = 0\}$ is traced. Thus, this method provides a design technique in identifying desirable values for the parameter $K$.

The root locus method allows one to identify the poles. As we studied in the relevant section, the pole information clearly identifies BIBO stability properties and all of the poles are to have negative real values for the closed-loop system to be BIBO stable.

Additionally, it lets one select desirable pole values: for example, poles with imaginary components lead to significant transient fluctuations and poles with real parts closer to the origin (on the left half plane) dominate the behaviour involving the response characteristics. A control engineer/designer may have reasons to choose certain poles over others.

For the approach to be applicable, the following key mathematical result is to be noted.

---

Nyquist Stability Criterion

With the Bode plot, we observed that we can identify the transfer function when a system is already stable. However, the Bode plot does not provide insights on how to adjust the controller so that stability is attained. The Root Locus method allows for parametrically adjusting the instability region. Complementing the Root Locus method, the Nyquist plot provides further insight on controller design and its robustness properties to parameter variations with also the phase of the transfer function appearing in the analysis (unlike the Root Locus method), to be discussed further below.

Recall first that a right-half plane pole of $1 + P(s)C(s)$ implies instability. In general, it is not difficult to identify the poles of the (open-loop) transfer function $P(s)C(s)$. Therefore, in the following we will assume that we know the number of right-half plane poles $P(s)C(s)$. Note also that the poles of $P(s)C(s)$ are the same as the poles of $1 + P(s)C(s)$; thus, we will assume that we know the number of right-half plane poles of $1 + P(s)C(s)$.

Constructing the Nyquist contour. For the Nyquist plot, we will construct a clockwise contour starting from $-iR$ to the origin and then to $+iR$ and then along a semi-circle of radius $R$ to close the curve. Later on we will take $R \to \infty$.

We will refer to this contour as a Nyquist contour.

If there is a pole of $L(s)$ on the imaginary axis, we should carefully exclude this from our contour: to exclude these, we divert the path along a semi-circle of a very small radius $r$ around the pole in a counter clock-wise fashion (which will later be taken to be arbitrarily close to 0). The exclusion of such a pole will not make a difference in the stability analysis: we focus on the zeroes of $1 + P(s)C(s)$ in the right-half plane (as such a pole cannot make $1 + P(s)C(s) = 0$).

Note: One can alternatively trace the contours counter-clockwise and then count $N$ as the number of counterclockwise encirclements. The result will be the same.

The proof builds on what is known as the principle of variation of the argument theorem, which we state and prove next.

Now, to apply this theorem, consider $f(s) = 1 + P(s)C(s)$, as noted earlier observe that the poles of $1 + P(s)C(s)$ are the same as the poles of $P(s)C(s)$. We are interested in the zeroes of $1 + P(s)C(s)$, and whether they are in the right-half plane. So, all we need to compute is the number of zeroes of $1 + P(s)C(s)$ through the difference given by the number of encirclements: Note now that the number of encirclements of $1 + P(s)C(s)$ around 0 is the same as the number of encirclements of $P(s)C(s)$ around $-1$. So, the number of zeroes in the right-half plane of $1 + P(s)C(s)$ will be $P$ plus the winding number. As a final note, in the Nyquist analysis, we construct the contour clockwise and apply the argument principle accordingly (so, the number of encirclements around $-1$ should be counted clock-wise). So, the number of interest is $N + P$, as claimed.

Computing the Nyquist plot. To compute the number of clock-wise encirclements, compute $L(i\omega)$ starting from $\omega = 0$ and increase $\omega$ to $\infty$. Observe that $L(i\omega) = \overline{L(-i\omega)}$, computing $L(i\omega)$ for $\omega > 0$ is sufficient to compute the values for $\omega < 0$. Finally, to compute $L(s)$ as $|s| \to \infty$, we note that often $|L(s)|$ converges to a constant as $|s| \to \infty$, and the essence of the encirclements is given by the changes observed as $s$ traces the imaginary axis.

Robustness via Nyquist plots

Nyquist's criterion suggests a robustness analysis, via gain and phase margins, mainly as a way to measure how far $P(s)L(s)$ is from $-1$ in both the magnitude (1) and phase ($\pi$) terms.

Roughly speaking, for systems which hit the real-line only once, in the complex plane the angle between $\pi$ and the location where the Nyquist plot hits the unit circle (i.e., with magnitude equaling 1) is called the phase stability margin.

The ratio between $-1$ and the point where the Nyquist plot hits the negative $x$-axis is called the gain stability margin.

Root locus method cannot give such a strong characterization for robustness. A useful implication is the small gain theorem presented in the following subsection.

System gain, passivity and the small gain theorem

Consider a linear system with feedback, which we assume to be stable. We generalize the observation above by viewing the input as one in $L_2(\mathbb{R}; \mathbb{C})$. Consider then the gain of a linear system with:

$\gamma := \sup_{u \in L_2(\mathbb{R}; \mathbb{C}): \|u\|_2 \neq 0} \frac{\|y\|_2}{\|u\|_2}$

We know, by Parseval's theorem (Theorem), that

$\gamma := \sup_{u \in L_2(\mathbb{R}; \mathbb{C}): \|u\|_2 \neq 0} \frac{\|\mathcal{F}_{CC}(y)\|_2}{\|\mathcal{F}_{CC}(u)\|_2}$

By writing $i\omega$ instead of $i2\pi f$ in the following, we have

$\mathcal{F}_{CC}(y)(\omega) = \mathcal{F}_{CC}(u)(\omega) G(i\omega),$

where $G$ is the closed-loop transfer function with $s = i\omega$. It can then be shown by noting that

$\int_\omega |\mathcal{F}_{CC}(u)(i\omega)|^2 |G(i\omega)|^2\, d\omega \le \left(\sup_\omega |G(i\omega)|^2\right) \int_\omega |\mathcal{F}_{CC}(u)(i\omega)|^2,$

the following holds:

$\gamma := \sup_{u \in L_2(\mathbb{R}; \mathbb{C}): \|u\|_2 \neq 0} \sqrt{\frac{\int_\omega |\mathcal{F}_{CC}(u)(i\omega) G(i\omega)|^2\, d\omega}{\|\mathcal{F}_{CC}(u)\|^2}} = \sup_{\omega \in \mathbb{R}} |G(i\omega)| =: \|G\|_\infty$

In the above, we write $\sup_{u \in L_2(\mathbb{R};\mathbb{C})}$, since a maximizing frequency $\omega^*$ exists, a corresponding input $u(t) = e^{i\omega^* t}$, will not be square integrable. However, this can be approximated arbitrarily well by truncation of the input and the output:

First observe that since $g$ is integrable, $G(i\omega)$ is continuous in $\omega$ (or $\hat{g}(2\pi f)$ is continuous in $f$). This can be shown by an application of the dominated convergence theorem. Let $\omega^* = 2\pi f^*$ and $u_K(t) = 1_{\{-\frac{K}{2} \le t < \frac{K}{2}\}} e^{i2\pi f^* t}$. As $K \to \infty$, the gain of the system will approximate $\gamma$ arbitrarily well. This is slightly beyond the scope of our course, but the reasoning is that the Fourier transform of $u_K(t) = 1_{\{-\frac{K}{2} \le t \le \frac{K}{2}\}} e^{i2\pi f^* t}$ will form an approximate identity sequence $\frac{\sin(\pi(f-f^*)K)}{\pi(f-f^*)}$, and $\frac{1}{K}|\hat{u}_K|^2(f)$ will also form an approximate identity-like sequence around $f^*$, since

$\frac{1}{K}\left(\frac{\sin(\pi(f-f^*))K)}{\pi(f-f^*)}\right)^2 = \frac{\sin(\pi(f-f^*))K)}{\pi(f-f^*)K} \frac{\sin(\pi(f-f^*))K)}{\pi(f-f^*)}$

Thus, the gain, when the input is $u_K$, is:

$\left(\frac{\int_f \big|\frac{\sin(\pi(f-f^*)K)}{\pi(f-f^*)}\big|^2 |\hat{g}(f)|^2\, df}{\int_f \big|\frac{\sin(\pi(f-f^*)K)}{\pi(f-f^*)}\big|^2\, df}\right)^{1/2}$

Now, by viewing the function

$R(f) := \frac{\big|\frac{\sin(\pi(f-f^*)K)}{\pi(f-f^*)}\big|^2}{\int_v \big|\frac{\sin(\pi(v-f^*)K)}{\pi(v-f^*)}\big|^2\, dv}$

as a probability density function which concentrates around $f^*$, we have that

$\left(\int_\omega \frac{\big|\frac{\sin(\pi(f-f^*)K)}{\pi(f-f^*)}\big|^2}{\int_\omega \big|\frac{\sin(\pi(f'-f^*)K)}{\pi(f'-f^*)}\big|^2\, df'} \hat{g}(f)^2\, df\right)$

converges to $|\hat{g}(f^*)|^2$ as $K$ gets larger and larger by the continuity of $\hat{g}$.

Thus, we can get arbitrarily close to the supremum gain by first approaching the supremum value over all frequencies arbitrarily well, and then by approaching the gain corresponding to such frequencies arbitrarily well with inputs constructed as above. Thus, we define the gain of a linear system as:

$\gamma := \|G(i\omega)\|_\infty = \sup_{\omega \in \mathbb{R}} |G(i\omega)| = \sup_{f \in \mathbb{R}} |\hat{g}(f)|$

Define a system to be $L_2$-stable if a bounded input, in the $L_2$-sense, leads to a bounded output in the $L_2$-sense. BIBO stability of a linear system implies $L_2$-stability: BIBO stability implies $\|h\|_1 < \infty$, which implies that $H(i\omega)$ is uniformly bounded for $\omega \in \mathbb{R}$ so that $|H(i\omega)| \le \|h\|_1$. Then, observe that the output $\|y^m\|_2^2 = \int |Y^m(i\omega)^2| \frac{1}{2\pi}\, d\omega = \int |H^2(i\omega)|^2 |U^m(i\omega)|^2 \frac{1}{2\pi}\, d\omega \le \|h\|_1^2 \int |(U^m(i\omega))^2| \frac{1}{2\pi}\, d\omega$. Thus, following Theorem, we can extend the domain (input function space) to be the entire $L_2$ space: Now, approximate any $u$ with $\|u\|_2 < \infty$ with the sequence $\{u^m \in \mathcal{S}, m \in \mathbb{N}\}$, taking the limit of the output sequence of which lets us conclude that the output $\|y\|_2$ is also bounded. Thus $\|y\|_2^2 = \int |Y(i\omega)^2| \frac{1}{2\pi}\, d\omega = \int |H^2(i\omega) U^2(i\omega)| \frac{1}{2\pi}\, d\omega < \infty$.

We next state a useful result on the verification of stability.

Note that $|H_1(i\omega)H_2(i\omega)| < 1$ uniformly for $\omega \in \mathbb{R}$ and therefore $\sup_\omega \frac{H_1(i\omega)}{1 + H_1(i\omega)H_2(i\omega)}$ is uniformly bounded from above by 1, since $1 + H_1(i\omega)H_2(i\omega)$ is uniformly bounded away from 0. The proof then follows from Nyquist's criterion: Since $H_1, H_2$ are stable, they have no poles in the right half plane. Furthermore, since $\gamma_1 \gamma_2 < 1$, $|C(s)P(s)|$ (there: $H_1(i\omega)H_2(i\omega)$) will be uniformly away from the point $-1$, thus a positive gain margin will be maintained and $-1$ will not be encircled. The closed-loop system is then stable.

The above concept does not involve phase properties, and it may be conservative for certain applications. On phase properties, there is a further commonly used concept of passivity: By Nyquist's criterion, for stable $P$ and $C$, if $P(s)C(s)$ is so that the phase is in $(-\pi, \pi)$ for all $s = i\omega, \omega \in \mathbb{R}$, then the closed loop system will be stable since $1 + P(s)C(s)$ will never hit the negative real axis and therefore $-1$ will not be encircled. Since $P$ and $C$ are stable, they have no poles in the right half-plane. The closed-loop system is then stable.

---

Predictive and Feedforward Control

In addition to stability and robustness, one may also wish to have desirable tracking and settling time properties (such as a low delay in the response). For such setups, it is often desirable to have a minimum-phase closed loop system (with no zeroes in the right half-plane). When one has a plant which is not minimum-phase, the design is often more challenging, due to the effective time delay present, as discussed in the relevant section.

In such cases, it may be necessary to go beyond classical PID control, and consider schemes such as predictive control or feedforward control where part of the control input can be computed by some pre-processing by possibly using feedback from the output.

To gain some insight, consider the following equation:

$\frac{P(s)C(s)}{1 + P(s)C(s)} = H_d(s)$

where $H_d(s)$ is a desired closed-loop system transfer function, which we may assume to be minimum-phase. Some algebra leads to:

$C(s) = \frac{1}{P(s)}\left(\frac{1}{H_d(s)} - 1\right)$

What we see is that if $P(s)$ has a zero in the right half-plane, the presence of $\frac{1}{P(s)}$ above leads to an unstable control transfer function.

For such settings, predictive control or feedforward control policies are often utilized in practice. Feedforward control is a general control paradigm where, unlike feedback control, the control is defined by an open-loop inner system in the control. The feedforward control is often used as a predictive correction term for the system. However, in many designs, the feedforward control is often added to a feedback control term leading to the overall control to the system with stronger robustness properties.

To gain further insight, let $P(s) = P_m(s)P_{nm}(s)$ where $P_m$ is minimum-phase and $P_{nm}$ is not minimum-phase. Consider $Y_f(s) = Y(s)P_{nm}^{-1}(s)$. This then effectively allows the control to only work with the minimum-phase portion of the system as if the feedback $Y_f$ only involves $P_m$. A PID control can then be designed for this system, and the output may be then subjected to $P_{nm}$ at the end: Observe that via $Y_f(s) = Y(s)P_{nm}^{-1}(s)$, $y_f$ acts as a prediction of the output $y$: For example, if $P_{nm}(s) = e^{-2s}$ were a pure delay system (see the Pade approximation discussion in the relevant section), then $Y_f(s) = Y(s)e^{2s}$ would essentially be written in time domain as $y_f(t) = y(t+2)$, serving as a prediction of the output. Thus, future is anticipated in the design.

For example, in reference tracking applications in which a path is given and the system is to closely steer the path, using feedback and causal reference information alone may lead to a weaker performance; here, future reference information can be used to compute an additional control input via an inverse analysis mapping reference paths to control actions. One can see that if it is possible to predict/estimate the future behaviour of a system to be tracked, this may make the application of a baseline control feasible (e.g. if an obstacle is known to be soon coming ahead, one may be able to slow down in anticipation). This additional control input is a predictive control and its presence can provide significant improvement in performance.

Some specific predictive control design methods are known as Smith Predictor Based or Internal Model Based Control Designs (which allow the control to anticipate the output and require further system model information), or Model Predictive Control (which is an optimization based control design and is more aligned with modern control theory).

---

The Routh-Hurwitz Stability Criterion

The root locus and Nyquist methods provide ways to analyze stability graphically. The Routh-Hurwitz criterion provides an algebraic test: given a polynomial (such as the characteristic polynomial of a closed-loop system), it determines whether all roots have negative real parts -- without computing the roots explicitly.

This is particularly useful for checking BIBO stability of a closed-loop transfer function $\frac{P(s)C(s)}{1 + P(s)C(s)}$, where the denominator polynomial $1 + P(s)C(s) = 0$ determines the pole locations.

Necessary Condition

Constructing the Routh Table

Given a polynomial $p(s) = a_n s^n + a_{n-1} s^{n-1} + \cdots + a_1 s + a_0$, the Routh table is constructed as follows. The first two rows are filled directly from the polynomial coefficients:

• Row $s^n$: $a_n, \; a_{n-2}, \; a_{n-4}, \; \ldots$
• Row $s^{n-1}$: $a_{n-1}, \; a_{n-3}, \; a_{n-5}, \; \ldots$
• Row $s^{n-2}$: $b_1, \; b_2, \; b_3, \; \ldots$
• Row $s^{n-3}$: $c_1, \; c_2, \; c_3, \; \ldots$
• Continue until row $s^0$.

Each subsequent entry is computed by a $2 \times 2$ determinant divided by the leading element of the row above:

$b_1 = \frac{a_{n-1} \cdot a_{n-2} - a_n \cdot a_{n-3}}{a_{n-1}}, \qquad b_2 = \frac{a_{n-1} \cdot a_{n-4} - a_n \cdot a_{n-5}}{a_{n-1}}$

and similarly for subsequent rows. The table has $n + 1$ rows total.

Worked Example

---

Exercises