add warning for nyquist_plot() when Nyquist criterion isn't met

murrayrm · murrayrm · commit 1e12a51f457a · 2022-04-15T17:52:12.000-07:00
diff --git a/control/freqplot.py b/control/freqplot.py
@@ -47,6 +47,7 @@
 import matplotlib.pyplot as plt
 import numpy as np
 import warnings
+from math import nan
 
 from .ctrlutil import unwrap
 from .bdalg import feedback
@@ -805,6 +806,28 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
         phase = -unwrap(np.angle(resp + 1))
         count = int(np.round(np.sum(np.diff(phase)) / np.pi, 0))
 
+        #
+        # Make sure that the enciriclements match the Nyquist criterion
+        #
+        # If the user specifies the frequency points to use, it is possible
+        # to miss enciriclements, so we check here to make sure that the
+        # Nyquist criterion is actually satisfied.
+        #
+        if isinstance(sys, (StateSpace, TransferFunction)):
+            P = (sys.pole().real > 0).sum() if indent_direction == 'right' \
+                else (sys.pole().real >= 0).sum()
+            Z = (sys.feedback().pole().real >= 0).sum()
+            if Z != count + P:
+                warnings.warn(
+                    "number of encirclements does not match Nyquist criterion;"
+                    " check frequency range and indent radius/direction",
+                    UserWarning, stacklevel=2)
+            elif indent_direction == 'none' and any(sys.pole().real == 0):
+                warnings.warn(
+                    "system has pure imaginary poles but indentation is"
+                    " turned off; results may be meaningless",
+                    RuntimeWarning, stacklevel=2)
+
         counts.append(count)
         contours.append(contour)
 
diff --git a/control/tests/nyquist_test.py b/control/tests/nyquist_test.py
@@ -98,9 +98,10 @@ def test_nyquist_basic():
     count, contour_indented = ct.nyquist_plot(
         sys, np.linspace(1e-4, 1e2, 100), return_contour=True)
     assert not all(contour_indented.real == 0)
-    count, contour = ct.nyquist_plot(
-        sys, np.linspace(1e-4, 1e2, 100), return_contour=True,
-        indent_direction='none')
+    with pytest.warns(UserWarning, match="encirclements does not match"):
+        count, contour = ct.nyquist_plot(
+            sys, np.linspace(1e-4, 1e2, 100), return_contour=True,
+            indent_direction='none')
     np.testing.assert_almost_equal(contour, 1j*np.linspace(1e-4, 1e2, 100))
 
     # Nyquist plot with poles at the origin, omega unspecified
@@ -166,10 +167,11 @@ def test_nyquist_fbs_examples():
 
     plt.figure()
     plt.title("Figure 10.10: L(s) = 3 (s+6)^2 / (s (s+1)^2) [zoom]")
-    count = ct.nyquist_plot(sys, omega_limits=[1.5, 1e3])
-    # Frequency limits for zoom give incorrect encirclement count
-    # assert _Z(sys) == count + _P(sys)
-    assert count == -1
+    with pytest.warns(UserWarning, match="encirclements does not match"):
+        count = ct.nyquist_plot(sys, omega_limits=[1.5, 1e3])
+        # Frequency limits for zoom give incorrect encirclement count
+        # assert _Z(sys) == count + _P(sys)
+        assert count == -1
 
 
 @pytest.mark.parametrize("arrows", [
@@ -276,8 +278,9 @@ def test_nyquist_indent_im():
 
     # Imaginary poles with no indentation
     plt.figure();
-    count = ct.nyquist_plot(
-        sys, np.linspace(0, 1e3, 1000), indent_direction='none')
+    with pytest.warns(UserWarning, match="encirclements does not match"):
+        count = ct.nyquist_plot(
+            sys, np.linspace(0, 1e3, 1000), indent_direction='none')
     plt.title(
         "Imaginary poles; indent_direction='none'; encirclements = %d" % count)
     assert _Z(sys) == count + _P(sys)
