Chapter 7: Fundamentals of Stability Analysis

Lesson 5: Special Cases — Repeated and Imaginary-Axis Roots

This lesson refines the stability concepts developed in previous lessons by focusing on repeated roots and imaginary-axis roots of the closed-loop characteristic polynomial. We analyze how root multiplicity and location on the imaginary axis affect time response, introduce special handling in the Routh–Hurwitz test, and show how these ideas appear in robotic joint models and their software implementations.

1. Root Multiplicity and Special Locations in the s-Plane

Consider a closed-loop characteristic polynomial \( p(s) \) of degree \( n \):

\[ p(s) = a_n s^n + a_{n-1} s^{n-1} + \dots + a_1 s + a_0, \quad a_n > 0. \]

Let the roots (poles) be \( s_1,\dots,s_r \) with multiplicities \( m_1,\dots,m_r \), so that

\[ p(s) = a_n \prod_{k=1}^{r} (s - s_k)^{m_k}, \quad \sum_{k=1}^{r} m_k = n. \]

We call a root repeated if \( m_k \ge 2 \). Two special geometric locations are:

  • Roots on the imaginary axis, \( s = \mathrm{j}\omega \) or \( s = -\mathrm{j}\omega \) with \( \omega \in \mathbb{R} \).
  • The origin, \( s = 0 \), which is part of the imaginary axis.

Stability classification (for continuous-time LTI systems) can be summarized as:

\[ \begin{aligned} &\text{Asymptotically stable} &\iff&\ \Re(s_k) < 0 \text{ for all } k. \\ &\text{Marginally (critically) stable} &\iff&\ \Re(s_k) \le 0 \text{ for all } k, \\ &&&\text{and all roots with } \Re(s_k) = 0 \text{ are simple (} m_k = 1\text{).} \\ &\text{Unstable} &\iff&\ \exists k \text{ with } \Re(s_k) > 0 \\ &&&\text{or some root on } \Re(s) = 0 \text{ has } m_k \ge 2. \end{aligned} \]

Thus, repeated roots on the imaginary axis are always a sign of instability, even though they are not in the open right half-plane.

flowchart TD
  P0["Characteristic polynomial p(s)"] --> P1["Compute roots s_k and multiplicities m_k"]
  P1 --> C1["Any root with Re(s_k) > 0 ?"]
  C1 -->|yes| U1["System unstable"]
  C1 -->|no| C2["Any root with Re(s_k) = 0 and m_k >= 2 ?"]
  C2 -->|yes| U2["System unstable \n(repeated imag-axis roots)"]
  C2 -->|no| C3["Any simple root with Re(s_k) = 0 ?"]
  C3 -->|yes| M1["Marginally stable \n(sustained oscillations)"]
  C3 -->|no| S1["Asymptotically stable"]
        

2. Time Response Induced by Repeated Poles

For an LTI system with transfer function \( G(s) = \dfrac{N(s)}{p(s)} \), inverse Laplace transform shows that each pole contributes a mode to the homogeneous solution. If \( s_0 \) is a pole of multiplicity \( m \), then partial fraction expansion yields terms of the form

\[ \frac{N(s)}{(s - s_0)^m} = \sum_{k=1}^{m} \frac{A_k}{(s - s_0)^k}, \quad A_k \in \mathbb{C}, \]

whose inverse Laplace transform is

\[ \mathcal{L}^{-1}\!\left\{\frac{A_k}{(s - s_0)^k}\right\} = \frac{A_k}{(k-1)!} t^{k-1} e^{s_0 t}. \]

Therefore the total contribution of the repeated pole \( s_0 \) is

\[ y_{s_0}(t) = \left( C_0 + C_1 t + \dots + C_{m-1} t^{m-1} \right) e^{s_0 t}, \quad C_i \in \mathbb{C}, \]

with stability behavior determined by \( \Re(s_0) \):

  • \( \Re(s_0) < 0 \): \( e^{s_0 t} \) decays exponentially, and the polynomial in \( t \) cannot prevent the mode from decaying. The system remains asymptotically stable (but with slower decay if \( m \) is large).
  • \( \Re(s_0) = 0 \): \( e^{s_0 t} \) has constant magnitude. If \( m = 1 \), the response is bounded (oscillatory). If \( m \ge 2 \), the polynomial in \( t \) causes unbounded growth.
  • \( \Re(s_0) > 0 \): Exponential growth dominates regardless of multiplicity; the system is unstable.

As an example, a robot joint axis with pure inertia \( J \) and no friction has plant model \( G(s) = \dfrac{1}{J s^2} \), which has a double pole at the origin. Even without right-half-plane poles, a step torque generates an unbounded ramp in angular velocity, illustrating instability due to repeated imaginary-axis roots.

3. Imaginary-Axis Roots and Marginal Stability

Now specialize to roots on the imaginary axis, \( s_0 = \mathrm{j}\omega_0 \) with \( \omega_0 \ne 0 \). By complex conjugate symmetry of real-coefficient polynomials, \( -\mathrm{j}\omega_0 \) is also a root with the same multiplicity.

For a simple pair of poles \( \pm \mathrm{j}\omega_0 \), the time response contains a term

\[ y_{\mathrm{osc}}(t) = C_1 \cos(\omega_0 t) + C_2 \sin(\omega_0 t), \]

which is bounded but does not decay. If all other poles satisfy \( \Re(s_k) < 0 \), the system is marginally stable.

For a repeated pair, e.g. poles of multiplicity 2 at \( \pm \mathrm{j}\omega_0 \), the response includes

\[ y_{\mathrm{osc}}(t) = \big( C_1 + C_2 t \big)\cos(\omega_0 t) + \big( C_3 + C_4 t \big)\sin(\omega_0 t), \]

whose envelope grows linearly with \( t \). Hence the system is unstable despite all poles lying on \( \Re(s) = 0 \). This subtle case is exactly why multiplicity matters in the stability definition.

From an engineering viewpoint, marginally stable systems with simple imaginary-axis poles are fragile: any small parameter variation (e.g. added actuator friction, flexibility, or controller gain change) can move these poles into the right half-plane or left half-plane, turning sustained oscillations into growing oscillations or decaying modes.

4. Routh–Hurwitz Criterion: Rows of Zeros and Imaginary-Axis Roots

In Lesson 3 we constructed the Routh array for \( p(s) \) and used sign changes in the first column to count right-half-plane roots. Special configurations of the Routh array indicate repeated and imaginary-axis roots.

4.1 Routh Array Recap

For brevity, let \( p(s) = a_4 s^4 + a_3 s^3 + a_2 s^2 + a_1 s + a_0 \) with \( a_4 > 0 \). The Routh array is

\[ \begin{array}{c|ccc} s^4 & a_4 & a_2 & a_0 \\ s^3 & a_3 & a_1 & 0 \\ s^2 & b_1 & b_2 & 0 \\ s^1 & c_1 & c_2 & 0 \\ s^0 & d_1 & 0 & 0 \end{array} \]

where, for example,

\[ b_1 = \frac{a_3 a_2 - a_4 a_1}{a_3}, \quad b_2 = \frac{a_3 a_0 - a_4 \cdot 0}{a_3} = \frac{a_3 a_0}{a_3} = a_0, \]

and similarly for the lower rows. All roots lie in the open left half-plane iff all entries in the first column are positive.

4.2 Zero First Element in a Row

If the first element of a row is zero (but the row is not entirely zero), e.g. \( c_1 = 0 \) while some other entry in the row is nonzero, the standard remedy is to replace \( c_1 \) by a small positive \( \varepsilon \), complete the array symbolically, and then let \( \varepsilon \to 0^+ \). The resulting sign changes in the first column reveal whether roots cross or lie on the imaginary axis.

This situation typically indicates a root on the imaginary axis and is often a warning that marginal stability or instability due to repeated imaginary-axis roots is present.

4.3 Entire Row of Zeros

If an entire row becomes zero, all elements in that row vanish. This signals that the polynomial has symmetrically located roots with respect to the imaginary axis, which may include imaginary-axis roots, possibly repeated.

Let the first all-zero row be associated with power \( s^k \). Then:

  1. Form the auxiliary polynomial \( A(s) \) from the row above (power \( s^{k+1} \)):

\[ A(s) = \alpha_0 s^{k+1} + \alpha_1 s^{k-1} + \alpha_2 s^{k-3} + \dots, \]

  1. Differentiate: \( A'(s) \) and use its coefficients to replace the zero row.
  2. Continue the Routh construction. The roots of \( A(s) = 0 \) include any imaginary-axis roots (possibly repeated).

Example. Consider \( p(s) = s^4 + 2 s^2 + 1 = (s^2 + 1)^2 \), which has double roots at \( s = \pm \mathrm{j} \). The Routh array starts with

\[ \begin{array}{c|ccc} s^4 & 1 & 2 & 1 \\ s^3 & 0 & 0 & 0 \\ \vdots & \vdots & \vdots & \vdots \end{array} \]

The \( s^3 \) row is entirely zero, so we form the auxiliary polynomial from the \( s^4 \) row: \( A(s) = s^4 + 2 s^2 + 1 \). Solving \( A(s) = 0 \) gives \( (s^2 + 1)^2 = 0 \), i.e. repeated imaginary-axis roots at \( s = \pm \mathrm{j} \). The associated time response contains terms growing like \( t \cos t \) and \( t \sin t \), confirming instability.

flowchart TD
  R0["Routh array construction"] --> R1["Any negative element in first column?"]
  R1 -->|yes| RU["Unstable (RHP roots)"]
  R1 -->|no| R2["Any row with first element 0 \nbut not all zeros?"]
  R2 -->|yes| RE["Replace first element \nby epsilon, analyze limit"]
  RE --> R3["Detect possible imag-axis roots"]
  R2 -->|no| R4["Any entire row of zeros?"]
  R4 -->|yes| RA["Form auxiliary polynomial A(s), \nuse A'(s), find imag-axis roots"]
  R4 -->|no| RS["All first-column entries > 0: \nasymptotically stable"]
        

5. Worked Example — Simple Imaginary-Axis Pair

Consider a fourth-order polynomial with a simple imaginary-axis pair and two stable real poles:

\[ p(s) = (s^2 + 1)(s + 1)(s + 2) = s^4 + 3 s^3 + 3 s^2 + 3 s + 2. \]

The roots are \( s = -1, -2, \pm \mathrm{j} \). The Routh array reveals a row of zeros at some intermediate step, leading (via the auxiliary polynomial) to the imaginary-axis roots \( \pm \mathrm{j} \). Since these roots are simple and all other roots lie strictly in the left half-plane, the system is marginally stable.

The homogeneous response contains:

  • Exponentially decaying terms from \( s = -1, -2 \).
  • A bounded oscillatory term \( C_1 \cos t + C_2 \sin t \) from \( \pm \mathrm{j} \).

Consequently, a robotic actuator described by such a characteristic polynomial would never settle to a constant steady state, but would continue to oscillate with constant amplitude. In practice, friction, actuator saturation, or additional damping usually shift these poles into the left half-plane and restore asymptotic stability (or into the right half-plane, causing growing oscillations, if the design is poor).

6. Python Lab — Detecting Repeated and Imaginary-Axis Roots

Python together with the numpy and python-control libraries is widely used in robotics for modeling servo loops and joint controllers. Below we demonstrate:

  • Computation of roots and multiplicities.
  • Classification into asymptotically stable, marginally stable, or unstable.
  • Application to a simple robot joint model.

import numpy as np
from control import tf, step_response  # python-control, common in robotics

def roots_with_multiplicity(den, tol=1e-6):
    """
    Compute roots of a real polynomial with approximate multiplicity.
    den: list of coefficients [a_n, ..., a_0]
    """
    r = np.roots(den)
    used = np.zeros(len(r), dtype=bool)
    clusters = []
    for i in range(len(r)):
        if used[i]:
            continue
        cluster = [r[i]]
        used[i] = True
        for j in range(i + 1, len(r)):
            if not used[j] and abs(r[j] - r[i]) < tol:
                used[j] = True
                cluster.append(r[j])
        clusters.append((np.mean(cluster), len(cluster)))
    return clusters  # list of (representative_root, multiplicity)

def classify_system(den, tol=1e-6):
    roots_mult = roots_with_multiplicity(den, tol)
    unstable = False
    marginal = False
    for s, m in roots_mult:
        if np.real(s) > tol:
            unstable = True
        elif abs(np.real(s)) <= tol:
            if m >= 2:
                unstable = True  # repeated imaginary-axis root
            else:
                marginal = True  # simple imag-axis root
    if unstable:
        return "unstable", roots_mult
    elif marginal:
        return "marginally stable", roots_mult
    else:
        return "asymptotically stable", roots_mult

# Example 1: (s^2 + 1)^2 => repeated imaginary-axis roots
den1 = [1, 0, 2, 0, 1]
cls1, info1 = classify_system(den1)
print("Example 1:", cls1)
for r, m in info1:
    print(" root =", r, " multiplicity =", m)

# Example 2: (s^2 + 1)(s + 1)(s + 2) => simple imag-axis pair
den2 = [1, 3, 3, 3, 2]
cls2, info2 = classify_system(den2)
print("\nExample 2:", cls2)
for r, m in info2:
    print(" root =", r, " multiplicity =", m)

# Robot joint with inertia J and viscous friction b
J = 0.01
b = 0.0   # set to zero => double integrator
num_joint = [1]
den_joint = [J, 0, 0]  # J s^2
cls3, info3 = classify_system(den_joint)
print("\nJoint example:", cls3)

# Build a transfer function and visualize step response
joint_sys = tf(num_joint, den_joint)
t, y = step_response(joint_sys)
print("Step response final value ~", y[-1])  # grows without bound for double integrator
      

In robotic applications, such classification can be integrated into automatic controller-tuning scripts to ensure that parameter updates (e.g., gains in a ROS-based joint controller) do not push the characteristic polynomial into unsafe regimes with repeated or right-half-plane roots.

7. C++ Lab — Companion Matrix and Stability Check (ROS-Friendly)

In C++, polynomial root-finding can be done via the eigenvalues of a companion matrix. This pattern integrates naturally with ROS and ROS 2 controllers that rely on the Eigen library for linear algebra.


#include <iostream>
#include <vector>
#include <complex>
#include <Eigen/Dense>

using namespace Eigen;

std::vector<std::pair<std::complex<double>, int> >
rootsWithMultiplicity(const std::vector<double> &den, double tol = 1e-6)
{
    int n = static_cast<int>(den.size()) - 1; // degree
    MatrixXd C = MatrixXd::Zero(n, n);

    // Companion matrix for a_n s^n + ... + a_0
    // First row: -a_{n-1}/a_n, ..., -a_0/a_n
    for (int j = 0; j < n; ++j)
        C(0, j) = -den[j + 1] / den[0]; // den[0] = a_n

    // Subdiagonal ones
    for (int i = 1; i < n; ++i)
        C(i, i - 1) = 1.0;

    EigenSolver<MatrixXd> es(C);
    VectorXcd r = es.eigenvalues();

    std::vector<bool> used(n, false);
    std::vector<std::pair<std::complex<double>, int> > clusters;

    for (int i = 0; i < n; ++i)
    {
        if (used[i]) continue;
        std::complex<double> ri = r(i);
        std::complex<double> sum = ri;
        int multiplicity = 1;
        used[i] = true;

        for (int j = i + 1; j < n; ++j)
        {
            if (!used[j] && std::abs(r(j) - ri) < tol)
            {
                used[j] = true;
                sum += r(j);
                ++multiplicity;
            }
        }
        clusters.emplace_back(sum / static_cast<double>(multiplicity),
                              multiplicity);
    }
    return clusters;
}

std::string classify(const std::vector<double> &den, double tol = 1e-6)
{
    auto roots = rootsWithMultiplicity(den, tol);
    bool unstable = false;
    bool marginal = false;
    for (auto &p : roots)
    {
        std::complex<double> s = p.first;
        int m = p.second;
        double re = s.real();
        if (re > tol)
            unstable = true;
        else if (std::abs(re) <= tol)
        {
            if (m >= 2) unstable = true;
            else marginal = true;
        }
    }
    if (unstable) return "unstable";
    if (marginal) return "marginally stable";
    return "asymptotically stable";
}

int main()
{
    // Example: (s^2 + 1)^2
    std::vector<double> den1 = {1.0, 0.0, 2.0, 0.0, 1.0};
    std::cout << "Example 1: " << classify(den1) << std::endl;

    // In a ROS controller, such a function could be called after
    // updating gains to verify that the closed-loop characteristic
    // polynomial remains acceptable.
    return 0;
}
      

When embedded inside a ROS node (e.g. alongside ros_control), such a check can be used to reject parameter updates that move the system into an unstable configuration characterized by repeated imaginary-axis roots or right-half-plane poles.

8. Java Lab — Apache Commons Math and Robotics Control Loops

Java is popular in educational and industrial robotics (e.g. FIRST Robotics with WPILib). Using Apache Commons Math, we can compute roots and classify stability for controller design.


import java.util.ArrayList;
import java.util.List;
import org.apache.commons.math3.analysis.polynomials.PolynomialFunction;
import org.apache.commons.math3.analysis.solvers.LaguerreSolver;
import org.apache.commons.math3.complex.Complex;

public class StabilityChecker {

    public static List<Complex[]> clusterRoots(Complex[] roots, double tol) {
        List<Complex[]> clusters = new ArrayList<>();
        boolean[] used = new boolean[roots.length];

        for (int i = 0; i < roots.length; ++i) {
            if (used[i]) continue;
            List<Complex> cluster = new ArrayList<>();
            cluster.add(roots[i]);
            used[i] = true;

            for (int j = i + 1; j < roots.length; ++j) {
                if (!used[j] && roots[j].subtract(roots[i]).abs() < tol) {
                    used[j] = true;
                    cluster.add(roots[j]);
                }
            }
            clusters.add(cluster.toArray(new Complex[0]));
        }
        return clusters;
    }

    public static String classify(double[] den, double tol) {
        PolynomialFunction p = new PolynomialFunction(den); // lowest order first
        LaguerreSolver solver = new LaguerreSolver();
        Complex[] roots = solver.solveAllComplex(den, 0.0);
        List<Complex[]> clusters = clusterRoots(roots, tol);

        boolean unstable = false;
        boolean marginal = false;

        for (Complex[] cluster : clusters) {
            Complex avg = Complex.ZERO;
            for (Complex c : cluster) {
                avg = avg.add(c);
            }
            avg = avg.divide(cluster.length);
            double re = avg.getReal();
            int mult = cluster.length;

            if (re > tol) unstable = true;
            else if (Math.abs(re) <= tol) {
                if (mult >= 2) unstable = true;
                else marginal = true;
            }
        }

        if (unstable) return "unstable";
        if (marginal) return "marginally stable";
        return "asymptotically stable";
    }

    public static void main(String[] args) {
        // Example: (s^2 + 1)^2 => coefficients in ascending powers of s
        double[] den1 = {1.0, 0.0, 2.0, 0.0, 1.0};
        System.out.println(classify(den1, 1e-6));
    }
}
      

In an FRC-style robot project, this stability checker can be integrated with Java-based control code (e.g., using WPILib) to monitor updated gain schedules and guard against introducing repeated imaginary-axis poles that could lead to unbounded oscillations.

9. MATLAB / Simulink and Wolfram Mathematica

9.1 MATLAB / Simulink

MATLAB with the Control System Toolbox and Robotics System Toolbox is a standard environment for robot manipulator control design. We can script detection of repeated and imaginary-axis roots and connect it to Simulink models of robot joints.


% Characteristic polynomial: (s^2 + 1)^2
den1 = [1 0 2 0 1];
r1 = roots(den1);

% Approximate multiplicities
tol = 1e-6;
clusters = {};
used = false(size(r1));
for i = 1:length(r1)
    if used(i), continue; end
    idx = find(abs(r1 - r1(i)) < tol);
    used(idx) = true;
    clusters{end+1} = r1(idx); %#ok<AGROW>
end

% Classify
unstable = false; marginal = false;
for k = 1:numel(clusters)
    c = clusters{k};
    s_avg = mean(c);
    m = numel(c);
    if real(s_avg) > tol
        unstable = true;
    elseif abs(real(s_avg)) <= tol
        if m >= 2, unstable = true; else marginal = true; end
    end
end

if unstable
    disp('unstable');
elseif marginal
    disp('marginally stable');
else
    disp('asymptotically stable');
end

% Simple robot joint model: inertia J, friction b
J = 0.01; b = 0.0;
num_joint = 1;
den_joint = [J 0 0];  % J s^2
G_joint = tf(num_joint, den_joint);
figure; step(G_joint); grid on;
title('Step response of double-integrator joint (unstable)');

% In Robotics System Toolbox or Peter Corke''s Robotics Toolbox,
% G_joint can serve as the linearized plant for a single revolute joint.
      

9.2 Wolfram Mathematica

Mathematica provides symbolic tools for exact stability analysis, which are helpful when deriving conditions for robot controller parameters symbolically.


(* Characteristic polynomial *)
p[s_] := (s^2 + 1)^2;

(* Roots and their multiplicities *)
roots = Roots[p[s] == 0, s];
(* roots is an algebraic representation; we can factor instead *)
Factor[p[s]]

(* Classification by real part and multiplicity *)
sol = Solve[p[s] == 0, s];
vals = s /. sol;

(* Count multiplicities *)
classes = Tally[vals];
classes // TableForm

(* Determine stability class *)
tol = 10^-8;
unstable = False; marginal = False;
Do[
  re = Re[c[[1]]];
  m  = c[[2]];
  If[re > tol, unstable = True,
    If[Abs[re] <= tol,
      If[m >= 2, unstable = True, marginal = True]
    ]
  ],
  {c, classes}
];

If[unstable, Print["unstable"],
  If[marginal, Print["marginally stable"],
    Print["asymptotically stable"]
  ]
];

(* In manipulator control design, p[s] can be obtained symbolically
   from the closed-loop transfer function of a linearized joint model,
   allowing exact parameter conditions for avoiding repeated imag-axis roots. *)
      

Mathematica can also interface with external robotics middleware (e.g. communicating with ROS nodes) to symbolically analyze linearized models of robot subsystems and then export numerically validated parameters that avoid repeated or right-half-plane poles.

10. Problems and Solutions

Problem 1 (Repeated Imaginary-Axis Poles): Consider the polynomial \( p(s) = s^4 + 2 s^2 + 1 \).
(a) Factor \( p(s) \).
(b) Classify the stability of a system with this characteristic polynomial.
(c) Sketch qualitatively the type of time response you expect for a bounded input.

Solution:

(a) We recognize a quadratic in \( s^2 \):

\[ p(s) = s^4 + 2 s^2 + 1 = (s^2 + 1)^2. \]

The roots are \( s = \pm \mathrm{j} \), each with multiplicity 2.

(b) All roots lie on the imaginary axis (no right-half-plane poles), but they are repeated. By the classification in Section 1, the system is unstable, because roots on \( \Re(s) = 0 \) with multiplicity \( m \ge 2 \) always produce unbounded growth in the response.

(c) The homogeneous response includes terms like \( t \cos t \) and \( t \sin t \), so the oscillation amplitude increases without bound. For a bounded input (e.g. a step), the total output is the sum of this unbounded homogeneous part and a bounded particular solution, hence still unbounded. Qualitatively, we expect oscillations whose envelope grows linearly in time.

Problem 2 (Robot Joint with Double Integrator): A simplified linear model of a robot joint axis with inertia \( J \) and negligible damping is \( G(s) = \dfrac{\theta(s)}{\tau(s)} = \dfrac{1}{J s^2} \), where \( \tau \) is torque and \( \theta \) is angle.
(a) Identify the poles and their multiplicity.
(b) Is the system asymptotically stable, marginally stable, or unstable?
(c) Interpret the result physically.

Solution:

(a) The denominator is \( J s^2 \), so the poles are \( s = 0 \) with multiplicity 2.

(b) The poles lie on the imaginary axis (\( \Re(s) = 0 \)) and are repeated (\( m = 2 \)). Therefore the system is unstable.

(c) Physically, applying a bounded torque (e.g. a step torque) produces a linearly increasing angular velocity and a quadratically increasing angle. The joint will spin faster and faster, which is unacceptable in most robotic applications. In practice, damping (mechanical or through feedback control) is added to move the poles into the left half-plane.

Problem 3 (Simple Imaginary-Axis Roots and Marginal Stability): Show that if the characteristic polynomial has a simple pair of purely imaginary roots \( \pm \mathrm{j}\omega_0 \) and all other roots satisfy \( \Re(s_k) < 0 \), then the system is marginally stable. Indicate what happens to the response as \( t \to \infty \).

Solution:

The contribution of \( \pm \mathrm{j}\omega_0 \) to the homogeneous solution is

\[ y_{\mathrm{osc}}(t) = C_1 \cos(\omega_0 t) + C_2 \sin(\omega_0 t), \]

which is bounded for all \( t \). All other modes decay exponentially because their poles lie in the open left half-plane. Hence for any initial condition and bounded input, the total response approaches a bounded periodic oscillation. The system does not converge to a fixed equilibrium, so it is not asymptotically stable; however, the response remains bounded, so it is not unstable either. This is precisely the definition of marginal stability.

Problem 4 (Routh Auxiliary Polynomial): For the polynomial \( p(s) = s^4 + 3 s^3 + 3 s^2 + 3 s + 2 \), show that the Routh array contains a row of zeros and find the auxiliary polynomial. Use it to identify the imaginary-axis roots.

Solution (outline):

  1. Form the first two rows: \( s^4: [1, 3, 2] \), \( s^3: [3, 3, 0] \).
  2. Using the standard Routh formula, construct the \( s^2 \) row. The next row (\( s^1 \)) turns out to be all zeros.
  3. The auxiliary polynomial is formed from the \( s^2 \) row: \( A(s) = s^4 + 2 s^2 + 1 \) (up to a positive scaling factor). Solving \( A(s) = 0 \) gives \( (s^2 + 1)^2 = 0 \), so the imaginary-axis roots are \( s = \pm \mathrm{j} \).
  4. Combining with the original factorization, \( p(s) = (s^2 + 1)(s + 1)(s + 2) \), we see that \( \pm \mathrm{j} \) are simple roots of \( p(s) \). Therefore the system is marginally stable.

Problem 5 (Design Guideline): Explain why in practice, control engineers usually avoid designs whose characteristic polynomials have simple imaginary-axis roots, even though such designs are marginally stable.

Solution:

Marginally stable systems are extremely sensitive to modeling errors, parameter drift, and nonlinearities such as saturation. Any small change in parameters can move the imaginary-axis poles into the right half-plane, leading to growing oscillations. In mechanical and robotic systems, nonlinear friction and backlash can also excite such modes. To guarantee robustness, engineers prefer designs that place all closed-loop poles strictly in the left half-plane with some margin (e.g. minimum real part), avoiding poles that lie exactly on the imaginary axis.

11. Summary

In this lesson we examined how root multiplicity and imaginary-axis location affect the stability of continuous-time LTI systems. Repeated poles generate polynomial factors in the time response; when these occur on the imaginary axis, the resulting oscillations grow without bound, leading to instability even in the absence of right-half-plane poles. Simple imaginary-axis poles yield bounded but non-decaying oscillations and correspond to marginally stable systems. We also saw how special patterns in the Routh array (zero first elements and entire rows of zeros) reveal these situations and how to resolve them via the auxiliary polynomial. Finally, we implemented computational checks in Python, C++, Java, MATLAB/Simulink, and Mathematica that are directly relevant to robotic control software, allowing automated detection and avoidance of these subtle but critical stability issues.

12. References

  1. Routh, E. J. (1877). A treatise on the stability of a given state of motion. Proceedings of the London Mathematical Society, 6, 86–97.
  2. Hurwitz, A. (1895). Über die Bedingungen, unter welchen eine Gleichung nur Wurzeln mit negativen reellen Teilen besitzt. Mathematische Annalen, 46(2), 273–284.
  3. Lyapunov, A. M. (1892). Stability of Motion. (English translation, Academic Press, 1966).
  4. Pontryagin, L. S. (1955). On the zeros of some elementary transcendental functions. AMS Translations, (2)1, 95–110.
  5. Jury, E. I. (1964). A stability criterion for linear discrete systems. Proceedings of the IRE, 50(6), 1493–1500.
  6. Tsypkin, Y. Z. (1964). Stability of motion in systems with delayed feedback. Automation and Remote Control, 25, 1–19.
  7. Zames, G. (1966). On the input-output stability of time-varying nonlinear feedback systems. Part I: Conditions derived using concepts of loop gain, conicity, and positivity. IEEE Transactions on Automatic Control, 11(2), 228–238.
  8. Bellman, R. E. (1949). On the stability of solutions of linear differential equations. Duke Mathematical Journal, 16(4), 627–647.
  9. Kalman, R. E. (1960). Contributions to the theory of optimal control. Boletín de la Sociedad Matemática Mexicana, 5, 102–119.
  10. Vinnicombe, G. (1993). Frequency domain uncertainty and the graph topology. IEEE Transactions on Automatic Control, 38(9), 1371–1383.