Consider a unity-feedback loop with open-loop transfer function \( L(s) = C(s)G(s) \), where \( G(s) \) is the plant and \( C(s) \) the controller. The closed-loop transfer function is

\[ T(s) = \frac{L(s)}{1 + L(s)}, \qquad \text{characteristic equation: } 1 + L(s) = 0. \]

Absolute stability asks whether all closed-loop poles satisfy \( \Re(s) < 0 \). In Nyquist terms (for systems satisfying the conditions of the Nyquist theorem), this is decided by the number of encirclements of the critical point \( -1 + j0 \).

Relative stability asks how far the system is from instability. For example:

• How much can the loop gain change before the system becomes unstable?
• How much additional phase lag (e.g., due to delays or filters) can we tolerate?
• How far are dominant closed-loop poles from the imaginary axis, affecting overshoot and settling time?

In the Nyquist diagram of \( L(j\omega) \), relative stability is primarily captured by how the Nyquist curve passes around and how close it approaches the point \( -1 + j0 \).

flowchart TD
  P["Plant model 'G(s)'"] --> C["Choose controller 'C(s)' and gain 'K'"]
  C --> L["Form open-loop 'L(s) = K C(s) G(s)'"]
  L --> NQ["Compute Nyquist curve 'L(j*omega)'"]
  NQ --> ENC["Check encirclements of '-1' (absolute stability)"]
  ENC --> RS["Measure distance to '-1' and margins (relative stability)"]
  RS --> D["Adjust 'K' or controller shape for desired performance"]
        

Throughout this lesson we assume that the Nyquist theorem conditions are satisfied and that the open-loop model has no poles on the imaginary axis (or that such poles are treated by standard indentation procedures covered previously).

A key geometric quantity in the Nyquist plane is the distance from the Nyquist curve to \( -1 \). For each frequency \( \omega \), define

\[ d(\omega) = \left| 1 + L(j\omega) \right| = \left| -1 - L(j\omega) \right|. \]

The minimum distance

\[ d_{\min} = \inf_{\omega \in \mathbb{R}} d(\omega) \]

measures how close the Nyquist curve comes to the critical point. Intuitively:

• Large \( d_{\min} \) → Nyquist curve comfortably far from \(-1\), strong relative stability.
• Small \( d_{\min} \) → Nyquist curve passes close to \(-1\), weak relative stability.

To connect this distance to robustness, consider a multiplicative perturbation of the loop:

\[ L_{\Delta}(s) = L(s)\,(1 + \Delta(s)), \qquad \|\Delta\|_{\infty} = \sup_{\omega}|\Delta(j\omega)|. \]

The perturbed closed-loop characteristic equation is

\[ 1 + L_{\Delta}(s) = 1 + L(s)(1 + \Delta(s)) = \underbrace{(1 + L(s))}_{\text{nominal}} + L(s)\Delta(s). \]

At a given frequency \( \omega \), writing \( z(\omega) = 1 + L(j\omega) \), the perturbation moves the point to \( z(\omega) + L(j\omega)\Delta(j\omega) \). A sufficient (but conservative) condition to avoid crossing the origin for all \( \omega \) is

\[ |L(j\omega)\Delta(j\omega)| < |1 + L(j\omega)| \quad \forall \omega \quad \Rightarrow \quad \| \Delta \|_{\infty} < \inf_{\omega} \frac{|1 + L(j\omega)|}{|L(j\omega)|}. \]

Since \( |1 + L(j\omega)| \) appears in the numerator, larger minimum distance \( d_{\min} \) guarantees robustness to larger modeling errors or gain variations.

flowchart TD
  NQ["Nyquist curve 'L(j*omega)'"] --> Pmin["Closest point 'z_min' to '-1'"]
  Pmin --> D["Compute 'd_min = |1 + L(j*omega)|'"]
  D --> ST["Large 'd_min' => strong relative stability"]
  D --> WT["Small 'd_min' => weak relative stability"]
        

Let the loop transfer function be decomposed as \( L(s) = K L_0(s) \), where \( K > 0 \) is an adjustable scalar gain (often the controller gain) and \( L_0(s) \) is fixed. The Nyquist plot scales radially with \( K \):

\[ L(j\omega) = K L_0(j\omega) \quad \Rightarrow \quad \text{for each } \omega,\; L(j\omega) \text{ lies on the ray from the origin through } L_0(j\omega). \]

As \( K \) increases from zero:

• The Nyquist curve initially shrinks to the origin (\( K \to 0 \)).
• For moderate \( K \), it moves outwards; at a critical gain, it may pass through \( -1 \).
• Beyond that critical gain, the number of encirclements of \( -1 \) changes and the closed loop may become unstable.

For a fixed frequency \( \omega^* \), the condition that the Nyquist curve passes through \( -1 \) is

\[ K L_0(j\omega^*) = -1. \]

Taking magnitudes and arguments, we obtain the two equations

\[ K |L_0(j\omega^*)| = 1, \qquad \arg L_0(j\omega^*) = -\pi + 2k\pi, \quad k \in \mathbb{Z}. \]

Any solution \( (\omega^*, K_{\text{crit}}) \) of these equations yields a critical gain \( K_{\text{crit}} \) at which a closed-loop pole lies on the imaginary axis. The Nyquist-based gain margin for nominal gain \( K_0 \) is then

\[ \text{GM} = \frac{K_{\text{crit}}}{K_0}, \]

provided that the Nyquist curve shows that the system is stable at \( K_0 \) and just unstable at \( K_{\text{crit}} \). A large gain margin indicates that the Nyquist curve is far from \( -1 \) in the radial direction, hence good relative stability against gain changes.

In Chapter 17 we will formalize gain (and phase) margins for general frequency responses, but here the Nyquist viewpoint already demonstrates how scaling \( K \) changes encirclements of \( -1 \).

For many servo and robotic joint systems, a dominant second-order closed-loop behavior is desired:

\[ s^2 + 2\zeta \omega_n s + \omega_n^2 = 0, \]

where \( \zeta \) is the damping ratio and \( \omega_n \) the natural frequency. One way to realize this closed-loop polynomial is to choose an open-loop transfer function

\[ L(s) = \frac{\omega_n^2}{s\,(s + 2\zeta \omega_n)}. \]

Indeed, the closed-loop characteristic equation \( 1 + L(s) = 0 \) becomes

\[ 1 + \frac{\omega_n^2}{s\,(s + 2\zeta \omega_n)} = 0 \quad \Longleftrightarrow \quad s^2 + 2\zeta \omega_n s + \omega_n^2 = 0. \]

Evaluating the Nyquist curve on the imaginary axis at \( s = j\omega_n \), we obtain

\[ L(j\omega_n) = \frac{\omega_n^2}{j\omega_n\left(j\omega_n + 2\zeta \omega_n\right)} = -\frac{1 + 2j\zeta}{4\zeta^2 + 1}. \]

Thus the Nyquist curve passes at the point

\[ \Re(L(j\omega_n)) = -\frac{1}{4\zeta^2 + 1}, \qquad \Im(L(j\omega_n)) = -\frac{2\zeta}{4\zeta^2 + 1}. \]

Both real and imaginary parts are negative, so the point lies in the third quadrant. As \( \zeta \) increases:

• The magnitude \( |L(j\omega_n)| = 1/\sqrt{4\zeta^2 + 1} \) decreases.
• The Nyquist curve is pulled towards the origin and away from \( -1 \).

Hence, larger damping ratio \( \zeta \) tends to increase the distance between the Nyquist curve and \( -1 \), improving relative stability and reducing overshoot. While the exact relationship between damping ratio and classical phase margin is somewhat involved, a key qualitative takeaway is:

• Nyquist curves that loop tightly around \( -1 \) correspond to small \( \zeta \) (underdamped response, high overshoot).
• Nyquist curves that stay away from \( -1 \) correspond to larger \( \zeta \) (well-damped response, modest overshoot).

Consider a third-order plant often used as a simplified model for a robotic joint with flexible dynamics:

\[ G(s) = \frac{1}{s(s+1)(s+2)}, \qquad C(s) = K. \]

The open-loop transfer function is \( L(s) = K G(s) \), and the closed-loop characteristic equation is

\[ 1 + L(s) = 0 \quad \Longleftrightarrow \quad s(s+1)(s+2) + K = 0 \quad \Longleftrightarrow \quad s^3 + 3s^2 + 2s + K = 0. \]

Using the Routh–Hurwitz criterion for a third-order polynomial \( s^3 + a_2 s^2 + a_1 s + a_0 \), asymptotic stability requires

\[ a_2 > 0,\quad a_1 > 0,\quad a_0 > 0,\quad a_2 a_1 - a_0 > 0. \]

For \( s^3 + 3s^2 + 2s + K \) we have \( a_2 = 3,\; a_1 = 2,\; a_0 = K \), so the conditions reduce to

\[ K > 0, \qquad 3\cdot 2 - K > 0 \quad \Rightarrow \quad 0 < K < 6. \]

Therefore:

• Absolute stability: the closed loop is stable for \( 0 < K < 6 \).
• Critical gain: \( K_{\text{crit}} = 6 \) produces one pole on the imaginary axis.

In the Nyquist diagram:

• For \( 0 < K < 6 \) the Nyquist curve does not encircle \( -1 \); the encirclement count agrees with Routh–Hurwitz.
• At \( K = 6 \), the Nyquist curve passes exactly through \( -1 \); \( d_{\min} = 0 \), and the system is marginally stable.
• For \( K > 6 \) the Nyquist curve encircles \( -1 \), and the closed loop is unstable.

Relative stability at a nominal gain \( K_0 \) can be gauged by how far the Nyquist curve for \( K_0 \) lies from the curve at \( K = 6 \). For example:

• \( K_0 = 1 \): the curve passes well away from \( -1 \), strong relative stability.
• \( K_0 = 4 \): the curve is closer to \( -1 \), weaker relative stability.
• \( K_0 = 5.5 \): the curve nearly touches \( -1 \), very weak relative stability.

In robotic joint control, a typical architecture uses nested loops: a fast inner current/torque loop and a slower outer position loop. For the outer loop, the effective plant often includes the inner loop dynamics, motor inertia, and gear train. The Nyquist plot of that effective plant with a proportional or PD position controller answers questions such as:

• How much can the stiffness (proportional gain) be increased to achieve tight tracking without inducing oscillations?
• How much additional delay (from computation, communication, or sensor filtering) can be tolerated?

In terms of Nyquist relative stability:

• A servo with a Nyquist curve that remains far from \( -1 \) is robust to changes in payload inertia and friction.
• A servo whose Nyquist curve nearly intersects \( -1 \) may become unstable when a heavier payload is attached or when extra filtering is added.

Thus, relative stability via Nyquist plots provides a geometric robustness check without explicitly recomputing closed-loop poles for each perturbation.

We now compute Nyquist plots and approximate relative stability metrics in Python using the python-control library. The same code structure can be embedded in robotic simulation frameworks such as roboticstoolbox-python when the plant is a linearized joint model.

import numpy as np
import matplotlib.pyplot as plt
import control as ct  # python-control

# Example plant: G(s) = 1 / (s (s + 1) (s + 2))
num = [1.0]
den = [1.0, 3.0, 2.0, 0.0]
G = ct.TransferFunction(num, den)

K = 2.0  # nominal proportional gain
L = K * G

# Nyquist plot
omega = np.logspace(-2, 2, 1000)
real, imag, freq = ct.nyquist_plot(L, omega, Plot=False)

plt.figure()
plt.plot(real, imag)
plt.plot(real, -imag)  # mirror for full Nyquist if needed
plt.scatter([-1.0], [0.0], marker="x")  # critical point -1
plt.axhline(0, linewidth=0.5)
plt.axvline(0, linewidth=0.5)
plt.xlabel("Re(L(j*omega))")
plt.ylabel("Im(L(j*omega))")
plt.title("Nyquist plot for L(s) with K = %.2f" % K)
plt.gca().set_aspect("equal", "box")
plt.grid(True)

# Approximate minimum distance to -1 for relative stability
z = real + 1j * imag
d = np.abs(1.0 + z)
d_min = np.min(d)
omega_min = freq[np.argmin(d)]
print("Approximate d_min =", d_min, "at omega =", omega_min)

# Approximate gain margin by scaling K until Nyquist passes -1
# (brute-force search over a grid of gains)
K_grid = np.linspace(0.1, 10.0, 200)
gm_est = None
for Kg in K_grid:
    Lg = Kg * G
    rg, ig, fg = ct.nyquist_plot(Lg, omega, Plot=False)
    zg = rg + 1j * ig
    dg = np.min(np.abs(1.0 + zg))
    if dg < 1e-2:  # "close enough" to -1
        gm_est = Kg / K
        break

print("Estimated gain margin (Nyquist-based) GM ~", gm_est)

plt.show()
      

In a robotic context, G can be the linearized joint dynamics from roboticstoolbox-python, and the above procedure provides an automated way to assess how far the Nyquist curve is from instability for different payloads or controller gains.

In C++, we often use Eigen for linear algebra and robotics libraries such as ROS-based stacks (ros_control) for integration. Here we sketch a simple Nyquist computation utility for a scalar transfer function, suitable for embedding into a robot joint controller analysis tool.

#include <iostream>
#include <complex>
#include <vector>
#include <cmath>

// Example: G(s) = 1 / (s (s + 1) (s + 2)), L(s) = K G(s)
std::complex<double> L_of_jw(double omega, double K) {
    std::complex<double> j(0.0, 1.0);
    std::complex<double> s = j * omega;
    std::complex<double> denom = s * (s + 1.0) * (s + 2.0);
    return K / denom;
}

int main() {
    double K = 2.0;
    std::vector<double> omega_vec;
    for (int k = 0; k <= 1000; ++k) {
        double omega = std::pow(10.0, -2.0 + 4.0 * k / 1000.0);
        omega_vec.push_back(omega);
    }

    double d_min = 1e9;
    double omega_min = 0.0;

    for (double omega : omega_vec) {
        std::complex<double> L = L_of_jw(omega, K);
        std::complex<double> z = 1.0 + L; // point in (1 + L(jw))-plane
        double d = std::abs(z);
        if (d < d_min) {
            d_min = d;
            omega_min = omega;
        }
    }

    std::cout << "Approx d_min = " << d_min
              << " at omega = " << omega_min << std::endl;

    // This data can be published to ROS topics or logged for plotting in Python/Matlab.
    return 0;
}
      

In robotic applications, L_of_jw can be extended to MIMO transfer matrices (using Eigen) and applied to each joint or to decoupled modes after linearization of a multi-DOF manipulator.

In Java, numerical utilities such as Apache Commons Math can be combined with robotics frameworks (e.g., FRC WPILib or custom simulation tools) to analyze Nyquist-based relative stability.

import org.apache.commons.math3.complex.Complex;

public class NyquistRelativeStability {

    // Example: G(s) = 1 / (s (s + 1) (s + 2)), L(s) = K G(s)
    public static Complex Lofjw(double omega, double K) {
        Complex j = new Complex(0.0, 1.0);
        Complex s = j.multiply(omega);
        Complex denom = s.multiply(s.add(1.0)).multiply(s.add(2.0));
        return new Complex(K, 0.0).divide(denom);
    }

    public static void main(String[] args) {
        double K = 2.0;
        double dMin = Double.POSITIVE_INFINITY;
        double omegaMin = 0.0;

        for (int k = 0; k <= 1000; ++k) {
            double omega = Math.pow(10.0, -2.0 + 4.0 * k / 1000.0);
            Complex L = Lofjw(omega, K);
            Complex z = L.add(1.0); // 1 + L(jw)
            double d = z.abs();
            if (d < dMin) {
                dMin = d;
                omegaMin = omega;
            }
        }

        System.out.println("Approx d_min = " + dMin + " at omega = " + omegaMin);

        // In a robotics codebase, this can be used to check stability of joint controllers
        // when gains are auto-tuned or adapted online.
    }
}
      

The same structure can be wrapped in a GUI or logging interface to explore how changing \( K \) moves the Nyquist curve relative to \( -1 \).

MATLAB with Control System Toolbox and Robotics System Toolbox is widely used in industrial and academic robotics. The following script computes Nyquist plots and relative stability measures.

% Example plant for a robotic joint axis
s = tf('s');
G = 1 / (s * (s + 1) * (s + 2));

K = 2;                 % nominal gain
L = K * G;

figure;
nyquist(L);
hold on;
plot(-1, 0, 'rx', 'MarkerSize', 8, 'LineWidth', 2);
grid on;
axis equal;
title(sprintf('Nyquist plot for L(s), K = %.2f', K));

% Compute gain and phase margins (Nyquist/Bode-consistent)
[GM, PM, Wcg, Wcp] = margin(L);
fprintf('Gain margin GM = %.2f (at w = %.3f rad/s)\n', GM, Wcg);
fprintf('Phase margin PM = %.2f deg (at w = %.3f rad/s)\n', PM, Wcp);

% Approximate d_min on a frequency grid
omega = logspace(-2, 2, 1000);
Lw = squeeze(freqresp(L, omega));
z = 1 + Lw;
d = abs(z);
[d_min, idx] = min(d);
fprintf('Approx d_min = %.4f at omega = %.3f rad/s\n', d_min, omega(idx));

% Simulink hint:
% - In Simulink, linearize a robot model using linmod or linearize (Robotics System Toolbox),
%   then use bode, nyquist, and margin on the linearized open-loop model.
      

For robotic manipulators, linearize can be applied around a given configuration to obtain \( G(s) \) for each joint, followed by Nyquist-based analysis of relative stability.

Mathematica provides symbolic and numerical tools for Nyquist plots and distance computations, which are useful when deriving analytic expressions for \( d_{\min} \) or exploring parameter dependence.

(* Define transfer function and loop gain *)
s = I*omega;
K = 2.0;
G[omega_] := 1.0 / (s*(s + 1)*(s + 2));
L[omega_] := K*G[omega];

(* Nyquist data *)
omegaRange = LogSpace[-2, 2, 1000]; (* define your own LogSpace or use Table *)
nyquistPoints = Table[L[w], {w, omegaRange}];

(* Distance to -1 *)
dList = Abs[1 + nyquistPoints];
dMin = Min[dList];
omegaMin = omegaRange[[First@First@Position[dList, dMin]]];

Print["Approx d_min = ", dMin, " at omega = ", omegaMin];

(* Nyquist plot *)
ListLinePlot[
  {ReIm /@ nyquistPoints},
  AxesLabel -> {"Re", "Im"},
  AspectRatio -> 1,
  PlotRange -> All,
  GridLines -> Automatic,
  Epilog -> {Red, PointSize[Large], Point[{-1, 0}]}
]
      

Symbolic manipulation (e.g., with Reduce or Solve) can be used to find analytic expressions for critical gains \( K_{\text{crit}} \) where the Nyquist curve intersects \( -1 \).

Problem 1 (Nyquist and Critical Gain): For the plant \( G(s) = 1/[s(s+1)(s+2)] \) with proportional controller \( C(s) = K \), show using Nyquist reasoning that the critical gain is \( K_{\text{crit}} = 6 \). Explain the connection with the Routh–Hurwitz result.

Solution: Nyquist stability requires that the number of encirclements \( N \) of \( -1 \) satisfies \( Z = N + P \), where \( P \) is the number of open-loop right-half-plane poles and \( Z \) the number of closed-loop right-half-plane poles. Here \( P = 0 \). For \( K \) very small, the Nyquist curve shrinks near the origin and clearly does not encircle \( -1 \), so \( N = 0 \) and the closed loop is stable. As \( K \) increases, the Nyquist curve scales radially. When it first passes through \( -1 \), the closed loop has a pole on the imaginary axis, i.e., marginal stability. That gain is \( K_{\text{crit}} \). Routh–Hurwitz applied to \( s^3 + 3s^2 + 2s + K = 0 \) gives the stability range \( 0 < K < 6 \), hence \( K_{\text{crit}} = 6 \). At this gain, the Nyquist curve passes through \( -1 \), consistent with the change in encirclement count as \( K \) crosses \( 6 \).

Problem 2 (Minimum Distance and Robust Gain Variations): Let \( L(s) = K_0 L_0(s) \) be stable in closed loop at nominal gain \( K_0 \). Suppose you measure \( d_{\min} = \inf_{\omega} |1 + L(j\omega)| > 0 \). Show that for small multiplicative gain perturbations \( L_{\Delta}(s) = L(s)(1+\delta) \) with constant \( \delta \in \mathbb{R} \), a sufficient condition for closed-loop stability is \( |\delta| < d_{\min}/\|L\|_{\infty} \), where \( \|L\|_{\infty} = \sup_{\omega}|L(j\omega)| \).

Solution: For the perturbed loop,

\[ 1 + L_{\Delta}(j\omega) = 1 + L(j\omega)(1+\delta) = (1 + L(j\omega)) + \delta L(j\omega). \]

If this quantity never vanishes for real \( \omega \), the Nyquist curve for \( L_{\Delta} \) never passes through the origin in the \( 1 + L_{\Delta} \)-plane and the encirclement count is unchanged, preserving stability. Using the reverse triangle inequality,

\[ |1 + L_{\Delta}(j\omega)| \ge |1 + L(j\omega)| - |\delta||L(j\omega)| \ge d_{\min} - |\delta| \|L\|_{\infty}. \]

A sufficient condition for \( |1 + L_{\Delta}(j\omega)| > 0 \) for all \( \omega \) is

\[ d_{\min} - |\delta| \|L\|_{\infty} > 0 \quad \Rightarrow \quad |\delta| < \frac{d_{\min}}{\|L\|_{\infty}}. \]

This shows that larger \( d_{\min} \) implies more robustness to gain perturbations.

Problem 3 (Qualitative Relationship Between Damping and Nyquist Shape): For the second-order closed-loop realization of Section 4, \( L(s) = \omega_n^2 / (s(s+2\zeta\omega_n)) \), qualitatively describe how the Nyquist curve changes as \( \zeta \) increases from small to large values. Discuss the implications for overshoot and relative stability.

Solution: As shown in Section 4, \( |L(j\omega_n)| = 1/\sqrt{4\zeta^2 + 1} \) decreases as \( \zeta \) increases. More generally, for each \( \omega \), the magnitude of \( L(j\omega) \) is reduced as \( \zeta \) grows. Therefore the Nyquist curve is pulled inward toward the origin, away from \( -1 \). For small \( \zeta \) (underdamped), the Nyquist curve extends further into the third quadrant and may approach \( -1 \) closely, indicating high overshoot and weak relative stability. For larger \( \zeta \), the curve is compact and far from \( -1 \), corresponding to reduced overshoot and stronger relative stability.

Problem 4 (Nyquist Check for a Simple Plant): Consider \( G(s) = 1/(s^2 + 2s + 2) \) with \( C(s) = K \). Show that the closed loop is stable for all \( K > 0 \), and explain how this appears in the Nyquist plot of \( L(s) = K G(s) \).

Solution: The characteristic equation is

\[ 1 + K G(s) = 0 \quad \Longleftrightarrow \quad s^2 + 2s + 2 + K = 0. \]

For \( s^2 + 2s + (2+K) \), the coefficients are positive for all \( K > 0 \), and the discriminant is \( \Delta = 4 - 4(2+K) = -4(1+K) < 0 \), so the poles are complex with real part \( -1 \). Thus \( \Re(s) = -1 \) for all \( K > 0 \), and the system is asymptotically stable. In the Nyquist diagram, \( L(s) = K G(s) \) scales radially with \( K \), but the Nyquist curve never passes through nor encircles \( -1 \) for positive \( K \). Relative stability with respect to gain changes is therefore very strong (infinite gain margin), although the damping (real part) remains fixed at \( -1 \).

Problem 5 (Design Hint from Nyquist Shape): A Nyquist plot for a unity-feedback servo loop passes near \( -1 \) and exhibits a large loop in the third quadrant. Suggest two controller modifications (without changing the plant) that would increase the distance from the Nyquist curve to \( -1 \), and briefly justify your answer.

Solution: Two possible modifications:

• Decrease loop gain in the frequency range where the Nyquist curve approaches \( -1 \) (e.g., reduce proportional gain or shape the loop with a lag compensator). This shrinks the Nyquist curve toward the origin, increasing its distance from \( -1 \).
• Add phase lead compensation around the gain-crossover frequency to provide positive phase shift and reduce the rotation of the Nyquist curve toward the negative real axis. This moves the curve away from the negative real axis near \( -1 \), again increasing the distance from \( -1 \). Both modifications enhance relative stability, at the cost of possible trade-offs in low-frequency tracking or high-frequency noise sensitivity.

In this lesson we extended the Nyquist stability criterion from a binary test of stability to a tool for assessing relative stability. We introduced the minimum distance \( d_{\min} = \inf_{\omega}|1 + L(j\omega)| \) and explained its role in robustness to gain and modeling uncertainties. We studied how gain scaling modifies encirclements of \( -1 \) and defined a Nyquist-based gain margin. For second-order closed-loop behavior we related damping ratio qualitatively to the shape of the Nyquist curve. Finally, we implemented Nyquist-based relative stability analysis across Python, C++, Java, MATLAB/Simulink, and Mathematica, with an eye toward robotic joint control and servo design.