From 18e997c037db361076c4fb9c61250ef7139d03f9 Mon Sep 17 00:00:00 2001
From: Richard Murray <murray@cds.caltech.edu>
Date: Sun, 10 Apr 2022 07:10:48 -0700
Subject: [PATCH 01/11] add warning for nyquist_plot() when indent_radius is
 too large

---
 control/freqplot.py           | 45 +++++++++++++++++++++++++++++++----
 control/tests/nyquist_test.py | 27 +++++++++++++++++++++
 2 files changed, 68 insertions(+), 4 deletions(-)

diff --git a/control/freqplot.py b/control/freqplot.py
index 7f29dce36..27a457d15 100644
--- a/control/freqplot.py
+++ b/control/freqplot.py
@@ -723,6 +723,7 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
                 and indent_direction != 'none':
             if sys.isctime():
                 splane_poles = sys.poles()
+                splane_cl_poles = sys.feedback().poles()
             else:
                 # map z-plane poles to s-plane, ignoring any at the origin
                 # because we don't need to indent for them
@@ -730,28 +731,64 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
                 zplane_poles = zplane_poles[~np.isclose(abs(zplane_poles), 0.)]
                 splane_poles = np.log(zplane_poles)/sys.dt
 
+                zplane_cl_poles = sys.feedback().poles()
+                zplane_cl_poles = zplane_cl_poles[
+                    ~np.isclose(abs(zplane_poles), 0.)]
+                splane_cl_poles = np.log(zplane_cl_poles)/sys.dt
+
+            #
+            # Check to make sure indent radius is small enough
+            #
+            # If there is a closed loop pole that is near the imaginary access
+            # at a point that is near an open loop pole, it is possible that
+            # indentation might skip or create an extraneous encirclement.
+            # We check for that situation here and generate a warning if that
+            # could happen.
+            #
+            for p_cl in splane_cl_poles:
+                # See if any closed loop poles are near the imaginary axis
+                if abs(p_cl.real) <= indent_radius:
+                    # See if any open loop poles are close to closed loop poles
+                    p_ol = splane_poles[
+                        (np.abs(splane_poles - p_cl)).argmin()]
+
+                    if abs(p_ol - p_cl) <= indent_radius:
+                        warnings.warn(
+                            "indented contour may miss closed loop pole; "
+                            "consider reducing indent_radius to be less than "
+                            f"{abs(p_ol - p_cl):5.2g}", stacklevel=2)
+
+            # See if we should add some frequency points near the origin
             if splane_contour[1].imag > indent_radius \
                     and np.any(np.isclose(abs(splane_poles), 0)) \
                     and not omega_range_given:
                 # add some points for quarter circle around poles at origin
+                # (these will get indented left or right below)
                 splane_contour = np.concatenate(
                     (1j * np.linspace(0., indent_radius, 50),
                     splane_contour[1:]))
+
             for i, s in enumerate(splane_contour):
                 # Find the nearest pole
                 p = splane_poles[(np.abs(splane_poles - s)).argmin()]
+
                 # See if we need to indent around it
                 if abs(s - p) < indent_radius:
+                    # Figure out how much to offset (simple trigonometry)
+                    offset = np.sqrt(indent_radius ** 2 - (s-p).imag ** 2) \
+                        -(s-p).real
+
+                    # Figure out which way to offset the contour point
                     if p.real < 0 or (np.isclose(p.real, 0) \
                             and indent_direction == 'right'):
                         # Indent to the right
-                        splane_contour[i] += \
-                            np.sqrt(indent_radius ** 2 - (s-p).imag ** 2)
+                        splane_contour[i] += offset
+
                     elif p.real > 0 or (np.isclose(p.real, 0) \
                             and indent_direction == 'left'):
                         # Indent to the left
-                        splane_contour[i] -= \
-                            np.sqrt(indent_radius ** 2 - (s-p).imag ** 2)
+                        splane_contour[i] -= offset
+
                     else:
                         ValueError("unknown value for indent_direction")
 
diff --git a/control/tests/nyquist_test.py b/control/tests/nyquist_test.py
index a001598a6..e93fef9ac 100644
--- a/control/tests/nyquist_test.py
+++ b/control/tests/nyquist_test.py
@@ -41,6 +41,33 @@ def test_nyquist_basic():
     N_sys = ct.nyquist_plot(sys)
     assert _Z(sys) == N_sys + _P(sys)
 
+    # Previously identified bug
+    #
+    # This example has an open loop pole at -0.06 and a closed loop pole at
+    # 0.06, so if you use an indent_radius of larger than 0.12, then the
+    # encirclements computed by nyquist_plot() will not properly predict
+    # stability.  A new warning messages was added to catch this case.
+    #
+    A = np.array([
+        [-3.56355873, -1.22980795, -1.5626527 , -0.4626829 , -0.16741484],
+        [-8.52361371, -3.60331459, -3.71574266, -0.43839201,  0.41893656],
+        [-2.50458726, -0.72361335, -1.77795489, -0.4038419 ,  0.52451147],
+        [-0.281183  ,  0.23391825,  0.19096003, -0.9771515 ,  0.66975606],
+        [-3.04982852, -1.1091943 , -1.40027242, -0.1974623 , -0.78930791]])
+    B = np.array([[-0.], [-1.42827213], [ 0.76806551], [-1.07987454], [0.]])
+    C = np.array([[-0.,  0.35557249,  0.35941791, -0., -1.42320969]])
+    D = np.array([[0]])
+    sys = ct.ss(A, B, C, D)
+
+    # With a small indent_radius, all should be fine
+    N_sys = ct.nyquist_plot(sys, indent_radius=0.001)
+    assert _Z(sys) == N_sys + _P(sys)
+
+    # With a larger indent_radius, we get a warning message + wrong answer
+    with pytest.warns(UserWarning, match="contour may miss closed loop pole"):
+        N_sys = ct.nyquist_plot(sys, indent_radius=0.2)
+        assert _Z(sys) != N_sys + _P(sys)
+
     # Unstable system
     sys = ct.tf([10], [1, 2, 2, 1])
     N_sys = ct.nyquist_plot(sys)

From 1e12a51f457a7ada030f5b80d25dabe27ad3f68a Mon Sep 17 00:00:00 2001
From: Richard Murray <murray@cds.caltech.edu>
Date: Sun, 10 Apr 2022 08:58:11 -0700
Subject: [PATCH 02/11] add warning for nyquist_plot() when Nyquist criterion
 isn't met

---
 control/freqplot.py           | 23 +++++++++++++++++++++++
 control/tests/nyquist_test.py | 21 ++++++++++++---------
 2 files changed, 35 insertions(+), 9 deletions(-)

diff --git a/control/freqplot.py b/control/freqplot.py
index 27a457d15..bfb2c7c6e 100644
--- 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
index e93fef9ac..fe4119666 100644
--- 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)

From 6da87066cf065bde32b23222519e653b65a34df7 Mon Sep 17 00:00:00 2001
From: Richard Murray <murray@cds.caltech.edu>
Date: Sun, 10 Apr 2022 14:32:12 -0700
Subject: [PATCH 03/11] add code to limit magnitude of response and enhance
 contour near poles * use plot_curve_magnitude to scale response at large
 magnitudes * add new line styles of scaled points on the curve * add points
 to contour near imaginary poles to avoid missing encirclements * add offsets
 to primary/mirror curve at large amplitude to avoid overlap * (still needs
 some unit tests, documention updates, and comparisons)

---
 control/freqplot.py           | 206 ++++++++++++++++++++++++++--------
 control/tests/nyquist_test.py |  60 +++++++---
 2 files changed, 200 insertions(+), 66 deletions(-)

diff --git a/control/freqplot.py b/control/freqplot.py
index bfb2c7c6e..d19a75ab1 100644
--- a/control/freqplot.py
+++ b/control/freqplot.py
@@ -519,11 +519,15 @@ def gen_zero_centered_series(val_min, val_max, period):
 
 # Default values for module parameter variables
 _nyquist_defaults = {
-    'nyquist.mirror_style': '--',
+    'nyquist.primary_style': ['-', ':'],        # style for primary curve
+    'nyquist.mirror_style': ['--', '-.'],       # style for mirror curve
     'nyquist.arrows': 2,
     'nyquist.arrow_size': 8,
-    'nyquist.indent_radius': 1e-1,
-    'nyquist.indent_direction': 'right',
+    'nyquist.indent_radius': 1e-6,              # indentation radius
+    'nyquist.indent_direction': 'right',        # indentation direction
+    'nyquist.indent_points': 50,                # number of points to insert
+    'nyquist.max_curve_magnitude': 20,
+    'nyquist.max_curve_offset': 0.02,           # percent offset of curves
 }
 
 
@@ -563,19 +567,32 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
     color : string
         Used to specify the color of the line and arrowhead.
 
-    mirror_style : string or False
-        Linestyle for mirror image of the Nyquist curve.  If `False` then
-        omit completely.  Default linestyle ('--') is determined by
-        config.defaults['nyquist.mirror_style'].
-
-    return_contour : bool
+    return_contour : bool, optional
         If 'True', return the contour used to evaluate the Nyquist plot.
 
-    label_freq : int
+    *args : :func:`matplotlib.pyplot.plot` positional properties, optional
+        Additional arguments for `matplotlib` plots (color, linestyle, etc)
+
+    **kwargs : :func:`matplotlib.pyplot.plot` keyword properties, optional
+        Additional keywords (passed to `matplotlib`)
+
+    Returns
+    -------
+    count : int (or list of int if len(syslist) > 1)
+        Number of encirclements of the point -1 by the Nyquist curve.  If
+        multiple systems are given, an array of counts is returned.
+
+    contour : ndarray (or list of ndarray if len(syslist) > 1)), optional
+        The contour used to create the primary Nyquist curve segment.  To
+        obtain the Nyquist curve values, evaluate system(s) along contour.
+
+    Additional Parameters
+    ---------------------
+    label_freq : int, optiona
         Label every nth frequency on the plot.  If not specified, no labels
         are generated.
 
-    arrows : int or 1D/2D array of floats
+    arrows : int or 1D/2D array of floats, optional
         Specify the number of arrows to plot on the Nyquist curve.  If an
         integer is passed. that number of equally spaced arrows will be
         plotted on each of the primary segment and the mirror image.  If a 1D
@@ -585,39 +602,51 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
         locations for the primary curve and the second row will be used for
         the mirror image.
 
-    arrow_size : float
+    arrow_size : float, optional
         Arrowhead width and length (in display coordinates).  Default value is
         8 and can be set using config.defaults['nyquist.arrow_size'].
 
-    arrow_style : matplotlib.patches.ArrowStyle
+    arrow_style : matplotlib.patches.ArrowStyle, optional
         Define style used for Nyquist curve arrows (overrides `arrow_size`).
 
-    indent_radius : float
-        Amount to indent the Nyquist contour around poles that are at or near
-        the imaginary axis.
-
-    indent_direction : str
+    indent_direction : str, optional
         For poles on the imaginary axis, set the direction of indentation to
         be 'right' (default), 'left', or 'none'.
 
-    warn_nyquist : bool, optional
-        If set to 'False', turn off warnings about frequencies above Nyquist.
-
-    *args : :func:`matplotlib.pyplot.plot` positional properties, optional
-        Additional arguments for `matplotlib` plots (color, linestyle, etc)
+    indent_points : int, optional
+        Number of points to insert in the Nyquist contour around poles that
+        are at or near the imaginary axis.
 
-    **kwargs : :func:`matplotlib.pyplot.plot` keyword properties, optional
-        Additional keywords (passed to `matplotlib`)
+    indent_radius : float, optional
+        Amount to indent the Nyquist contour around poles that are at or near
+        the imaginary axis.
 
-    Returns
-    -------
-    count : int (or list of int if len(syslist) > 1)
-        Number of encirclements of the point -1 by the Nyquist curve.  If
-        multiple systems are given, an array of counts is returned.
+    max_curve_magnitude : float, optional
+        Restrict the maximum magnitude of the Nyquist plot to this value.
+        Portions of the Nyquist plot whose magnitude is restricted are
+        plotted using a different line style.
+
+    max_curve_offset : float, optional
+        When plotting scaled portion of the Nyquist plot, increase/decrease
+        the magnitude by this fraction of the max_curve_magnitude to allow
+        any overlaps between the primary and mirror curves to be avoided.
+
+    mirror_style : [str, str] or False
+        Linestyles for mirror image of the Nyquist curve.  The first element
+        is used for unscaled portions of the Nyquist curve, the second element
+        is used for portions that are scaled (using max_curve_magnitude).  If
+        `False` then omit completely.  Default linestyle (['--', '-.']) is
+        determined by config.defaults['nyquist.mirror_style'].
+
+    primary_style : [str, str]
+        Linestyles for primary image of the Nyquist curve.  The first element
+        is used for unscaled portions of the Nyquist curve, the second
+        element is used for scaled portions that are scaled (using
+        max_curve_magnitude). Default linestyle (['-', ':']) is determined by
+        config.defaults['nyquist.mirror_style'].
 
-    contour : ndarray (or list of ndarray if len(syslist) > 1)), optional
-        The contour used to create the primary Nyquist curve segment.  To
-        obtain the Nyquist curve values, evaluate system(s) along contour.
+    warn_nyquist : bool, optional
+        If set to 'False', turn off warnings about frequencies above Nyquist.
 
     Notes
     -----
@@ -668,8 +697,6 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
 
     # Get values for params (and pop from list to allow keyword use in plot)
     omega_num = config._get_param('freqplot', 'number_of_samples', omega_num)
-    mirror_style = config._get_param(
-        'nyquist', 'mirror_style', kwargs, _nyquist_defaults, pop=True)
     arrows = config._get_param(
         'nyquist', 'arrows', kwargs, _nyquist_defaults, pop=True)
     arrow_size = config._get_param(
@@ -679,6 +706,28 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
         'nyquist', 'indent_radius', kwargs, _nyquist_defaults, pop=True)
     indent_direction = config._get_param(
         'nyquist', 'indent_direction', kwargs, _nyquist_defaults, pop=True)
+    indent_points = config._get_param(
+        'nyquist', 'indent_points', kwargs, _nyquist_defaults, pop=True)
+    max_curve_magnitude = config._get_param(
+        'nyquist', 'max_curve_magnitude', kwargs, _nyquist_defaults, pop=True)
+    max_curve_offset = config._get_param(
+        'nyquist', 'max_curve_offset', kwargs, _nyquist_defaults, pop=True)
+
+    # Set line styles for the curves
+    def _parse_linestyle(style_name, allow_false=False):
+        style = config._get_param(
+            'nyquist', style_name, kwargs, _nyquist_defaults, pop=True)
+        if isinstance(style, str):
+            # Only one style provided, use the default for the other
+            style = [style, _nyquist_defaults['nyquist.' + style_name][1]]
+        if (allow_false and style is False) or \
+           (isinstance(style, list) and len(style) == 2):
+            return style
+        else:
+            raise ValueError(f"invalid '{style_name}': {style}")
+
+    primary_style = _parse_linestyle('primary_style')
+    mirror_style = _parse_linestyle('mirror_style', allow_false=True)
 
     # If argument was a singleton, turn it into a tuple
     if not isinstance(syslist, (list, tuple)):
@@ -759,15 +808,39 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
                             "consider reducing indent_radius to be less than "
                             f"{abs(p_ol - p_cl):5.2g}", stacklevel=2)
 
-            # See if we should add some frequency points near the origin
-            if splane_contour[1].imag > indent_radius \
-                    and np.any(np.isclose(abs(splane_poles), 0)) \
-                    and not omega_range_given:
-                # add some points for quarter circle around poles at origin
+            #
+            # See if we should add some frequency points near imaginary poles
+            #
+            for p in splane_poles:
+                # See if we need to process this pole (skip any that is on
+                # the not near or on the negative omega axis + user override)
+                if p.imag < 0 or abs(p.real) > indent_radius or \
+                   omega_range_given:
+                    continue
+
+                # Find the frequencies before the pole frequency
+                below_points = np.argwhere(
+                    splane_contour.imag - abs(p.imag) < -indent_radius)
+                if below_points.size > 0:
+                    first_point = below_points[-1].item()
+                    start_freq = p.imag - indent_radius
+                else:
+                    # Add the points starting at the beginning of the contour
+                    assert splane_contour[0] == 0
+                    first_point = 0
+                    start_freq = 0
+
+                above_points = np.argwhere(
+                    splane_contour.imag - abs(p.imag) > indent_radius)
+                last_point = above_points[0].item()
+
+                # Add points for half/quarter circle around pole frequency
                 # (these will get indented left or right below)
-                splane_contour = np.concatenate(
-                    (1j * np.linspace(0., indent_radius, 50),
-                    splane_contour[1:]))
+                splane_contour = np.concatenate((
+                    splane_contour[0:first_point+1],
+                    (1j * np.linspace(
+                        start_freq, p.imag + indent_radius, indent_points)),
+                    splane_contour[last_point:]))
 
             for i, s in enumerate(splane_contour):
                 # Find the nearest pole
@@ -849,19 +922,55 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
                 arrow_style = mpl.patches.ArrowStyle(
                     'simple', head_width=arrow_size, head_length=arrow_size)
 
-            # Save the components of the response
-            x, y = resp.real, resp.imag
+            # Find the different portions of the curve (with scaled pts marked)
+            reg_mask = abs(resp) > max_curve_magnitude
+            scale_mask = ~reg_mask \
+                & np.concatenate((~reg_mask[1:], ~reg_mask[-1:])) \
+                & np.concatenate((~reg_mask[0:1], ~reg_mask[:-1]))
+
+            # Rescale the points with large magnitude
+            resp[reg_mask] /= (np.abs(resp[reg_mask]) / max_curve_magnitude)
 
-            # Plot the primary curve
-            p = plt.plot(x, y, '-', color=color, *args, **kwargs)
+            # Plot the regular portions of the curve (and grab the color)
+            x_reg = np.ma.masked_where(reg_mask, resp.real)
+            y_reg = np.ma.masked_where(reg_mask, resp.imag)
+            p = plt.plot(
+                x_reg, y_reg, primary_style[0], color=color, *args, **kwargs)
             c = p[0].get_color()
+
+            # Plot the scaled sections of the curve (changing linestyle)
+            x_scl = np.ma.masked_where(scale_mask, resp.real)
+            y_scl = np.ma.masked_where(scale_mask, resp.imag)
+            plt.plot(
+                x_scl * (1 + max_curve_offset), y_scl * (1 + max_curve_offset),
+                primary_style[1], color=c, *args, **kwargs)
+
+            # Plot the primary curve (invisible) for setting arrows
+            x, y = resp.real.copy(), resp.imag.copy()
+            x[reg_mask] *= (1 + max_curve_offset)
+            y[reg_mask] *= (1 + max_curve_offset)
+            p = plt.plot(x, y, linestyle='None', color=c, *args, **kwargs)
+
+            # Add arrows
             ax = plt.gca()
             _add_arrows_to_line2D(
                 ax, p[0], arrow_pos, arrowstyle=arrow_style, dir=1)
 
             # Plot the mirror image
             if mirror_style is not False:
-                p = plt.plot(x, -y, mirror_style, color=c, *args, **kwargs)
+                # Plot the regular and scaled segments
+                plt.plot(
+                    x_reg, -y_reg, mirror_style[0], color=c, *args, **kwargs)
+                plt.plot(
+                    x_scl * (1 - max_curve_offset),
+                    -y_scl * (1 - max_curve_offset),
+                    mirror_style[1], color=c, *args, **kwargs)
+
+                # Add the arrows (on top of an invisible contour)
+                x, y = resp.real.copy(), resp.imag.copy()
+                x[reg_mask] *= (1 - max_curve_offset)
+                y[reg_mask] *= (1 - max_curve_offset)
+                p = plt.plot(x, -y, linestyle='None', color=c, *args, **kwargs)
                 _add_arrows_to_line2D(
                     ax, p[0], arrow_pos, arrowstyle=arrow_style, dir=-1)
 
@@ -982,7 +1091,6 @@ def _add_arrows_to_line2D(
 #
 # Gang of Four plot
 #
-
 # TODO: think about how (and whether) to handle lists of systems
 def gangof4_plot(P, C, omega=None, **kwargs):
     """Plot the "Gang of 4" transfer functions for a system
diff --git a/control/tests/nyquist_test.py b/control/tests/nyquist_test.py
index fe4119666..4f142e901 100644
--- a/control/tests/nyquist_test.py
+++ b/control/tests/nyquist_test.py
@@ -94,14 +94,18 @@ def test_nyquist_basic():
         contour[contour.real == 0],
         1j*np.linspace(0, 1e2, 100)[contour.real == 0])
 
+    #
     # Make sure that we can turn off frequency modification
+    #
+    # Start with a case where indentation should occur
     count, contour_indented = ct.nyquist_plot(
-        sys, np.linspace(1e-4, 1e2, 100), return_contour=True)
+        sys, np.linspace(1e-4, 1e2, 100), indent_radius=1e-2,
+        return_contour=True)
     assert not all(contour_indented.real == 0)
     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')
+            sys, np.linspace(1e-4, 1e2, 100), indent_radius=1e-2,
+            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
@@ -115,15 +119,20 @@ def test_nyquist_basic():
     assert _Z(sys) == count + _P(sys)
 
     # Nyquist plot with poles on imaginary axis, omega specified
+    # (can miss encirclements due to the imaginary poles at +/- 1j)
     sys = ct.tf([1], [1, 3, 2]) * ct.tf([1], [1, 0, 1])
-    count = ct.nyquist_plot(sys, np.linspace(1e-3, 1e1, 1000))
-    assert _Z(sys) == count + _P(sys)
+    with pytest.warns(UserWarning, match="does not match") as records:
+        count = ct.nyquist_plot(sys, np.linspace(1e-3, 1e1, 1000))
+    if len(records) == 0:
+        assert _Z(sys) == count + _P(sys)
 
     # Nyquist plot with poles on imaginary axis, omega specified, with contour
     sys = ct.tf([1], [1, 3, 2]) * ct.tf([1], [1, 0, 1])
-    count, contour = ct.nyquist_plot(
-        sys, np.linspace(1e-3, 1e1, 1000), return_contour=True)
-    assert _Z(sys) == count + _P(sys)
+    with pytest.warns(UserWarning, match="does not match") as records:
+        count, contour = ct.nyquist_plot(
+            sys, np.linspace(1e-3, 1e1, 1000), return_contour=True)
+    if len(records) == 0:
+        assert _Z(sys) == count + _P(sys)
 
     # Nyquist plot with poles on imaginary axis, return contour
     sys = ct.tf([1], [1, 3, 2]) * ct.tf([1], [1, 0, 1])
@@ -229,10 +238,10 @@ def test_nyquist_indent_default(indentsys):
 def test_nyquist_indent_dont(indentsys):
     # first value of default omega vector was 0.1, replaced by 0. for contour
     # indent_radius is larger than 0.1 -> no extra quater circle around origin
-    count, contour = ct.nyquist_plot(indentsys,
-                                     plot=False,
-                                     indent_radius=.1007,
-                                     return_contour=True)
+    with pytest.warns(UserWarning, match="encirclements does not match"):
+        count, contour = ct.nyquist_plot(
+            indentsys, omega=[0, 0.2, 0.3, 0.4], indent_radius=.1007,
+            plot=False, return_contour=True)
     np.testing.assert_allclose(contour[0], .1007+0.j)
     # second value of omega_vector is larger than indent_radius: not indented
     assert np.all(contour.real[2:] == 0.)
@@ -240,9 +249,8 @@ def test_nyquist_indent_dont(indentsys):
 
 def test_nyquist_indent_do(indentsys):
     plt.figure();
-    count, contour = ct.nyquist_plot(indentsys,
-                                     indent_radius=0.01,
-                                     return_contour=True)
+    count, contour = ct.nyquist_plot(
+        indentsys, indent_radius=0.01, return_contour=True)
     plt.title("Pole at origin; indent_radius=0.01; encirclements = %d" % count)
     assert _Z(indentsys) == count + _P(indentsys)
     # indent radius is smaller than the start of the default omega vector
@@ -305,6 +313,19 @@ def test_nyquist_exceptions():
     with pytest.warns(UserWarning, match="above Nyquist"):
         ct.nyquist_plot(sys, np.logspace(-2, 3))
 
+def test_linestyle_checks():
+    sys = ct.rss(2, 1, 1)
+
+    # Things that should work
+    ct.nyquist_plot(sys, primary_style=['-', '-'], mirror_style=['-', '-'])
+    ct.nyquist_plot(sys, mirror_style=None)
+
+    with pytest.raises(ValueError, match="invalid 'primary_style'"):
+        ct.nyquist_plot(sys, primary_style=False)
+
+    with pytest.raises(ValueError, match="invalid 'mirror_style'"):
+        ct.nyquist_plot(sys, mirror_style=0.2)
+
 
 if __name__ == "__main__":
     #
@@ -333,11 +354,16 @@ def test_nyquist_exceptions():
     test_nyquist_encirclements()
 
     print("Indentation checks")
-    test_nyquist_indent()
+    s = ct.tf('s')
+    indentsys = 3 * (s+6)**2 / (s * (s+1)**2)
+    test_nyquist_indent_default(indentsys)
+    test_nyquist_indent_do(indentsys)
+    test_nyquist_indent_left(indentsys)
 
     print("Unusual Nyquist plot")
     sys = ct.tf([1], [1, 3, 2]) * ct.tf([1], [1, 0, 1])
     plt.figure()
-    plt.title("Poles: %s" % np.array2string(sys.poles(), precision=2, separator=','))
+    plt.title("Poles: %s" %
+              np.array2string(sys.poles(), precision=2, separator=','))
     count = ct.nyquist_plot(sys)
     assert _Z(sys) == count + _P(sys)

From 3d3a5593742c938f58a69a11a64406eadcd5dfd5 Mon Sep 17 00:00:00 2001
From: Richard Murray <murray@cds.caltech.edu>
Date: Sun, 10 Apr 2022 23:01:22 -0700
Subject: [PATCH 04/11] smooth curve offsets to avoid discontinuities

---
 control/freqplot.py | 90 +++++++++++++++++++++++++++++++++++++++------
 1 file changed, 79 insertions(+), 11 deletions(-)

diff --git a/control/freqplot.py b/control/freqplot.py
index d19a75ab1..05ed7cf4f 100644
--- a/control/freqplot.py
+++ b/control/freqplot.py
@@ -519,8 +519,8 @@ def gen_zero_centered_series(val_min, val_max, period):
 
 # Default values for module parameter variables
 _nyquist_defaults = {
-    'nyquist.primary_style': ['-', ':'],        # style for primary curve
-    'nyquist.mirror_style': ['--', '-.'],       # style for mirror curve
+    'nyquist.primary_style': ['-', '-.'],        # style for primary curve
+    'nyquist.mirror_style': ['--', ':'],       # style for mirror curve
     'nyquist.arrows': 2,
     'nyquist.arrow_size': 8,
     'nyquist.indent_radius': 1e-6,              # indentation radius
@@ -696,6 +696,7 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
         kwargs.pop('arrow_length', False)
 
     # Get values for params (and pop from list to allow keyword use in plot)
+    omega_num_given = omega_num is not None
     omega_num = config._get_param('freqplot', 'number_of_samples', omega_num)
     arrows = config._get_param(
         'nyquist', 'arrows', kwargs, _nyquist_defaults, pop=True)
@@ -736,8 +737,13 @@ def _parse_linestyle(style_name, allow_false=False):
     omega, omega_range_given = _determine_omega_vector(
         syslist, omega, omega_limits, omega_num)
     if not omega_range_given:
-        # Start contour at zero frequency
-        omega[0] = 0.
+        if omega_num_given:
+            # Just reset the starting point
+            omega[0] = 0.0
+        else:
+            # Insert points between the origin and the first frequency point
+            omega = np.concatenate((
+                np.linspace(0, omega[0], indent_points), omega[1:]))
 
     # Go through each system and keep track of the results
     counts, contours = [], []
@@ -938,17 +944,23 @@ def _parse_linestyle(style_name, allow_false=False):
                 x_reg, y_reg, primary_style[0], color=color, *args, **kwargs)
             c = p[0].get_color()
 
+            # Figure out how much to offset the curve: the offset goes from
+            # zero at the start of the scaled section to max_curve_offset as
+            # we move along the curve
+            curve_offset = _compute_curve_offset(
+                resp, scale_mask, max_curve_offset)
+
             # Plot the scaled sections of the curve (changing linestyle)
             x_scl = np.ma.masked_where(scale_mask, resp.real)
             y_scl = np.ma.masked_where(scale_mask, resp.imag)
             plt.plot(
-                x_scl * (1 + max_curve_offset), y_scl * (1 + max_curve_offset),
+                x_scl * (1 + curve_offset), y_scl * (1 + curve_offset),
                 primary_style[1], color=c, *args, **kwargs)
 
             # Plot the primary curve (invisible) for setting arrows
             x, y = resp.real.copy(), resp.imag.copy()
-            x[reg_mask] *= (1 + max_curve_offset)
-            y[reg_mask] *= (1 + max_curve_offset)
+            x[reg_mask] *= (1 + curve_offset[reg_mask])
+            y[reg_mask] *= (1 + curve_offset[reg_mask])
             p = plt.plot(x, y, linestyle='None', color=c, *args, **kwargs)
 
             # Add arrows
@@ -962,14 +974,14 @@ def _parse_linestyle(style_name, allow_false=False):
                 plt.plot(
                     x_reg, -y_reg, mirror_style[0], color=c, *args, **kwargs)
                 plt.plot(
-                    x_scl * (1 - max_curve_offset),
-                    -y_scl * (1 - max_curve_offset),
+                    x_scl * (1 - curve_offset),
+                    -y_scl * (1 - curve_offset),
                     mirror_style[1], color=c, *args, **kwargs)
 
                 # Add the arrows (on top of an invisible contour)
                 x, y = resp.real.copy(), resp.imag.copy()
-                x[reg_mask] *= (1 - max_curve_offset)
-                y[reg_mask] *= (1 - max_curve_offset)
+                x[reg_mask] *= (1 - curve_offset[reg_mask])
+                y[reg_mask] *= (1 - curve_offset[reg_mask])
                 p = plt.plot(x, -y, linestyle='None', color=c, *args, **kwargs)
                 _add_arrows_to_line2D(
                     ax, p[0], arrow_pos, arrowstyle=arrow_style, dir=-1)
@@ -1087,6 +1099,62 @@ def _add_arrows_to_line2D(
         arrows.append(p)
     return arrows
 
+#
+# Function to compute Nyquist curve offsets
+#
+# This function computes a smoothly varying offset that starts and ends at
+# zero at the ends of a scaled segment.
+#
+def _compute_curve_offset(resp, mask, max_offset):
+    # Compute the arc length along the curve
+    s_curve = np.cumsum(
+        np.sqrt(np.diff(resp.real) ** 2 + np.diff(resp.imag) ** 2))
+
+    # Initialize the offset
+    offset = np.zeros(resp.size)
+    arclen = np.zeros(resp.size)
+
+    # Walk through the response and keep track of each continous component
+    i, nsegs = 0, 0
+    while i < resp.size:
+        # Skip the regular segment
+        while i < resp.size and mask[i]:
+            i += 1              # Increment the counter
+            if i == resp.size:
+                break;
+            # Keep track of the arclength
+            arclen[i] = arclen[i-1] + np.abs(resp[i] - resp[i-1])
+
+        nsegs += 0.5
+        if i == resp.size:
+            break;
+
+        # Save the starting offset of this segment
+        seg_start = i
+
+        # Walk through the scaled segment
+        while i < resp.size and not mask[i]:
+            i += 1
+            if i == resp.size:  # See if we are done with this segment
+                break;
+            # Keep track of the arclength
+            arclen[i] = arclen[i-1] + np.abs(resp[i] - resp[i-1])
+
+        nsegs += 0.5
+        if i == resp.size:
+            break;
+
+        # Save the ending offset of this segment
+        seg_end = i
+
+        # Now compute the scaling for this segment
+        s_segment = arclen[seg_end-1] - arclen[seg_start]
+        offset[seg_start:seg_end] = max_offset * s_segment/s_curve[-1] * \
+            np.sin(np.pi * (arclen[seg_start:seg_end]
+                            - arclen[seg_start])/s_segment)
+
+    return offset
+
 
 #
 # Gang of Four plot

From aef3b0618faa2af26372219ea0c20e0ea34a16b4 Mon Sep 17 00:00:00 2001
From: Richard Murray <murray@cds.caltech.edu>
Date: Tue, 12 Apr 2022 23:04:22 -0700
Subject: [PATCH 05/11] tweak max_curve_magnitude + add start point

---
 control/freqplot.py           | 54 +++++++++++++++++++++++++++--------
 control/tests/nyquist_test.py |  7 +++++
 2 files changed, 49 insertions(+), 12 deletions(-)

diff --git a/control/freqplot.py b/control/freqplot.py
index 05ed7cf4f..3b16b935c 100644
--- a/control/freqplot.py
+++ b/control/freqplot.py
@@ -519,15 +519,17 @@ def gen_zero_centered_series(val_min, val_max, period):
 
 # Default values for module parameter variables
 _nyquist_defaults = {
-    'nyquist.primary_style': ['-', '-.'],        # style for primary curve
-    'nyquist.mirror_style': ['--', ':'],       # style for mirror curve
-    'nyquist.arrows': 2,
-    'nyquist.arrow_size': 8,
-    'nyquist.indent_radius': 1e-6,              # indentation radius
+    'nyquist.primary_style': ['-', '-.'],       # style for primary curve
+    'nyquist.mirror_style': ['--', ':'],        # style for mirror curve
+    'nyquist.arrows': 2,                        # number of arrors around curve
+    'nyquist.arrow_size': 8,                    # pixel size for arrows
+    'nyquist.indent_radius': 1e-4,              # indentation radius
     'nyquist.indent_direction': 'right',        # indentation direction
     'nyquist.indent_points': 50,                # number of points to insert
-    'nyquist.max_curve_magnitude': 20,
-    'nyquist.max_curve_offset': 0.02,           # percent offset of curves
+    'nyquist.max_curve_magnitude': 20,          # clip large values
+    'nyquist.max_curve_offset': 0.02,           # offset of primary/mirror
+    'nyquist.start_marker': 'o',                # marker at start of curve
+    'nyquist.start_marker_size': 4,             # size of the maker
 }
 
 
@@ -618,8 +620,9 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
         are at or near the imaginary axis.
 
     indent_radius : float, optional
-        Amount to indent the Nyquist contour around poles that are at or near
-        the imaginary axis.
+        Amount to indent the Nyquist contour around poles on or near the
+        imaginary axis. Portions of the Nyquist plot corresponding to indented
+        portions of the contour are plotted using a different line style.
 
     max_curve_magnitude : float, optional
         Restrict the maximum magnitude of the Nyquist plot to this value.
@@ -638,13 +641,22 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
         `False` then omit completely.  Default linestyle (['--', '-.']) is
         determined by config.defaults['nyquist.mirror_style'].
 
-    primary_style : [str, str]
+    primary_style : [str, str], optional
         Linestyles for primary image of the Nyquist curve.  The first element
         is used for unscaled portions of the Nyquist curve, the second
         element is used for scaled portions that are scaled (using
         max_curve_magnitude). Default linestyle (['-', ':']) is determined by
         config.defaults['nyquist.mirror_style'].
 
+    start_marker : str, optional
+        Matplotlib marker to use to mark the starting point of the Nyquist
+        plot.  Defaults value is 'o' and can be set using
+        config.defaults['nyquist.start_marker'].
+
+    start_marker_size : float, optional
+        Start marker size (in display coordinates).  Default value is
+        4 and can be set using config.defaults['nyquist.start_marker_size'].
+
     warn_nyquist : bool, optional
         If set to 'False', turn off warnings about frequencies above Nyquist.
 
@@ -713,6 +725,10 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
         'nyquist', 'max_curve_magnitude', kwargs, _nyquist_defaults, pop=True)
     max_curve_offset = config._get_param(
         'nyquist', 'max_curve_offset', kwargs, _nyquist_defaults, pop=True)
+    start_marker = config._get_param(
+        'nyquist', 'start_marker', kwargs, _nyquist_defaults, pop=True)
+    start_marker_size = config._get_param(
+        'nyquist', 'start_marker_size', kwargs, _nyquist_defaults, pop=True)
 
     # Set line styles for the curves
     def _parse_linestyle(style_name, allow_false=False):
@@ -848,6 +864,7 @@ def _parse_linestyle(style_name, allow_false=False):
                         start_freq, p.imag + indent_radius, indent_points)),
                     splane_contour[last_point:]))
 
+            # Indent points that are too close to a pole
             for i, s in enumerate(splane_contour):
                 # Find the nearest pole
                 p = splane_poles[(np.abs(splane_poles - s)).argmin()]
@@ -929,13 +946,21 @@ def _parse_linestyle(style_name, allow_false=False):
                     'simple', head_width=arrow_size, head_length=arrow_size)
 
             # Find the different portions of the curve (with scaled pts marked)
-            reg_mask = abs(resp) > max_curve_magnitude
+            reg_mask = np.logical_or(
+                np.abs(resp) > max_curve_magnitude,
+                contour.real != 0)
+            # reg_mask = np.logical_or(
+            #     np.abs(resp.real) > max_curve_magnitude,
+            #     np.abs(resp.imag) > max_curve_magnitude)
+
             scale_mask = ~reg_mask \
                 & np.concatenate((~reg_mask[1:], ~reg_mask[-1:])) \
                 & np.concatenate((~reg_mask[0:1], ~reg_mask[:-1]))
 
             # Rescale the points with large magnitude
-            resp[reg_mask] /= (np.abs(resp[reg_mask]) / max_curve_magnitude)
+            rescale = np.logical_and(
+                reg_mask, abs(resp) > max_curve_magnitude)
+            resp[rescale] *= max_curve_magnitude / abs(resp[rescale])
 
             # Plot the regular portions of the curve (and grab the color)
             x_reg = np.ma.masked_where(reg_mask, resp.real)
@@ -986,6 +1011,11 @@ def _parse_linestyle(style_name, allow_false=False):
                 _add_arrows_to_line2D(
                     ax, p[0], arrow_pos, arrowstyle=arrow_style, dir=-1)
 
+            # Mark the start of the curve
+            if start_marker:
+                plt.plot(resp[0].real, resp[0].imag, start_marker,
+                         color=c, markersize=start_marker_size)
+
             # Mark the -1 point
             plt.plot([-1], [0], 'r+')
 
diff --git a/control/tests/nyquist_test.py b/control/tests/nyquist_test.py
index 4f142e901..730b36696 100644
--- a/control/tests/nyquist_test.py
+++ b/control/tests/nyquist_test.py
@@ -360,6 +360,13 @@ def test_linestyle_checks():
     test_nyquist_indent_do(indentsys)
     test_nyquist_indent_left(indentsys)
 
+    # Generate a figuring showing effects of different parameters
+    sys = 3 * (s+6)**2 / (s * (s**2 + 1e-4 * s + 1))
+    plt.figure()
+    ct.nyquist_plot(sys)
+    ct.nyquist_plot(sys, max_curve_magnitude=15)
+    ct.nyquist_plot(sys, indent_radius=1e-6, max_curve_magnitude=25)
+
     print("Unusual Nyquist plot")
     sys = ct.tf([1], [1, 3, 2]) * ct.tf([1], [1, 0, 1])
     plt.figure()

From c25c0cb39a00c22078d83a288302a21e4fc9e634 Mon Sep 17 00:00:00 2001
From: Richard Murray <murray@cds.caltech.edu>
Date: Wed, 13 Apr 2022 22:38:07 -0700
Subject: [PATCH 06/11] add warning if encirclements is not near an integer

---
 control/freqplot.py            | 24 +++++++++++++++++++++++-
 control/tests/descfcn_test.py  |  4 ++--
 control/tests/freqresp_test.py |  5 +++--
 control/tests/nyquist_test.py  | 11 +++++++++++
 4 files changed, 39 insertions(+), 5 deletions(-)

diff --git a/control/freqplot.py b/control/freqplot.py
index 3b16b935c..7710c0d53 100644
--- a/control/freqplot.py
+++ b/control/freqplot.py
@@ -523,6 +523,7 @@ def gen_zero_centered_series(val_min, val_max, period):
     'nyquist.mirror_style': ['--', ':'],        # style for mirror curve
     'nyquist.arrows': 2,                        # number of arrors around curve
     'nyquist.arrow_size': 8,                    # pixel size for arrows
+    'nyquist.encirclement_threshold': 0.05,     # warning threshold
     'nyquist.indent_radius': 1e-4,              # indentation radius
     'nyquist.indent_direction': 'right',        # indentation direction
     'nyquist.indent_points': 50,                # number of points to insert
@@ -611,6 +612,11 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
     arrow_style : matplotlib.patches.ArrowStyle, optional
         Define style used for Nyquist curve arrows (overrides `arrow_size`).
 
+    encirclement_threshold : float, optional
+        Define the threshold for generating a warning if the number of net
+        encirclements is a non-integer value.  Default value is 0.05 and can
+        be set using config.defaults['nyquist.encirclement_threshold'].
+
     indent_direction : str, optional
         For poles on the imaginary axis, set the direction of indentation to
         be 'right' (default), 'left', or 'none'.
@@ -717,6 +723,9 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
     arrow_style = config._get_param('nyquist', 'arrow_style', kwargs, None)
     indent_radius = config._get_param(
         'nyquist', 'indent_radius', kwargs, _nyquist_defaults, pop=True)
+    encirclement_threshold = config._get_param(
+        'nyquist', 'encirclement_threshold', kwargs,
+        _nyquist_defaults, pop=True)
     indent_direction = config._get_param(
         'nyquist', 'indent_direction', kwargs, _nyquist_defaults, pop=True)
     indent_points = config._get_param(
@@ -750,8 +759,11 @@ def _parse_linestyle(style_name, allow_false=False):
     if not isinstance(syslist, (list, tuple)):
         syslist = (syslist,)
 
+    # Determine the range of frequencies to use, based on args/features
     omega, omega_range_given = _determine_omega_vector(
         syslist, omega, omega_limits, omega_num)
+
+    # If omega was not specified explicitly, start at omega = 0
     if not omega_range_given:
         if omega_num_given:
             # Just reset the starting point
@@ -852,6 +864,7 @@ def _parse_linestyle(style_name, allow_false=False):
                     first_point = 0
                     start_freq = 0
 
+                # Find the frequencies after the pole frequency
                 above_points = np.argwhere(
                     splane_contour.imag - abs(p.imag) > indent_radius)
                 last_point = above_points[0].item()
@@ -900,7 +913,15 @@ def _parse_linestyle(style_name, allow_false=False):
 
         # Compute CW encirclements of -1 by integrating the (unwrapped) angle
         phase = -unwrap(np.angle(resp + 1))
-        count = int(np.round(np.sum(np.diff(phase)) / np.pi, 0))
+        encirclements = np.sum(np.diff(phase)) / np.pi
+        count = int(np.round(encirclements, 0))
+
+        # Let the user know if the count might not make sense
+        if abs(encirclements - count) > encirclement_threshold:
+            warnings.warn(
+                "number of encirclements was a non-integer value; this can"
+                " happen is contour is not closed, possibly based on a"
+                " frequency range that does not include zero.")
 
         #
         # Make sure that the enciriclements match the Nyquist criterion
@@ -1534,6 +1555,7 @@ def _determine_omega_vector(syslist, omega_in, omega_limits, omega_num,
                                     num=omega_num, endpoint=True)
     else:
         omega_out = np.copy(omega_in)
+
     return omega_out, omega_range_given
 
 
diff --git a/control/tests/descfcn_test.py b/control/tests/descfcn_test.py
index 796ad9034..0ebfda446 100644
--- a/control/tests/descfcn_test.py
+++ b/control/tests/descfcn_test.py
@@ -140,7 +140,7 @@ def test_describing_function(fcn, amin, amax):
 def test_describing_function_plot():
     # Simple linear system with at most 1 intersection
     H_simple = ct.tf([1], [1, 2, 2, 1])
-    omega = np.logspace(-1, 2, 100)
+    omega = np.logspace(-2, 2, 100)
 
     # Saturation nonlinearity
     F_saturation = ct.descfcn.saturation_nonlinearity(1)
@@ -160,7 +160,7 @@ def test_describing_function_plot():
 
     # Multiple intersections
     H_multiple = H_simple * ct.tf(*ct.pade(5, 4)) * 4
-    omega = np.logspace(-1, 3, 50)
+    omega = np.logspace(-2, 3, 50)
     F_backlash = ct.descfcn.friction_backlash_nonlinearity(1)
     amp = np.linspace(0.6, 5, 50)
     xsects = ct.describing_function_plot(H_multiple, F_backlash, amp, omega)
diff --git a/control/tests/freqresp_test.py b/control/tests/freqresp_test.py
index 0e35a38ea..573fd6359 100644
--- a/control/tests/freqresp_test.py
+++ b/control/tests/freqresp_test.py
@@ -81,8 +81,9 @@ def test_nyquist_basic(ss_siso):
         tf_siso, plot=False, return_contour=True, omega_num=20)
     assert len(contour) == 20
 
-    count, contour = nyquist_plot(
-        tf_siso, plot=False, omega_limits=(1, 100), return_contour=True)
+    with pytest.warns(UserWarning, match="encirclements was a non-integer"):
+        count, contour = nyquist_plot(
+            tf_siso, plot=False, omega_limits=(1, 100), return_contour=True)
     assert_allclose(contour[0], 1j)
     assert_allclose(contour[-1], 100j)
 
diff --git a/control/tests/nyquist_test.py b/control/tests/nyquist_test.py
index 730b36696..3b383f655 100644
--- a/control/tests/nyquist_test.py
+++ b/control/tests/nyquist_test.py
@@ -219,6 +219,17 @@ def test_nyquist_encirclements():
     plt.title("Pole at the origin; encirclements = %d" % count)
     assert _Z(sys) == count + _P(sys)
 
+    # Non-integer number of encirclements
+    plt.figure();
+    sys = 1 / (s**2 + s + 1)
+    with pytest.warns(UserWarning, match="encirclements was a non-integer"):
+        count = ct.nyquist_plot(sys, omega_limits=[0.5, 1e3])
+    with pytest.warns(None) as records:
+        count = ct.nyquist_plot(
+            sys, omega_limits=[0.5, 1e3], encirclement_threshold=0.2)
+    assert len(records) == 0
+    plt.title("Non-integer number of encirclements [%g]" % count)
+
 
 @pytest.fixture
 def indentsys():

From 67c1f3f68944feecae6e9028744f2c007746423b Mon Sep 17 00:00:00 2001
From: Richard Murray <murray@cds.caltech.edu>
Date: Thu, 14 Apr 2022 18:59:33 -0700
Subject: [PATCH 07/11] add legacy settings for nyquist

---
 control/config.py             |  9 +++++++++
 control/tests/nyquist_test.py | 11 +++++++++++
 2 files changed, 20 insertions(+)

diff --git a/control/config.py b/control/config.py
index 605fbcb23..32f5f2eef 100644
--- a/control/config.py
+++ b/control/config.py
@@ -267,6 +267,15 @@ def use_legacy_defaults(version):
     #
     reset_defaults()            # start from a clean slate
 
+    # Version 0.9.2:
+    if major == 0 and minor < 9 or (minor == 9 and patch < 2):
+        from math import inf
+
+        # Reset Nyquist defaults
+        set_defaults('nyquist', indent_radius=0.1, max_curve_magnitude=inf,
+                     max_curve_offset=0, primary_style=['-', '-'],
+                     mirror_style=['--', '--'], start_marker_size=0)
+
     # Version 0.9.0:
     if major == 0 and minor < 9:
         # switched to 'array' as default for state space objects
diff --git a/control/tests/nyquist_test.py b/control/tests/nyquist_test.py
index 3b383f655..8ac083535 100644
--- a/control/tests/nyquist_test.py
+++ b/control/tests/nyquist_test.py
@@ -337,6 +337,17 @@ def test_linestyle_checks():
     with pytest.raises(ValueError, match="invalid 'mirror_style'"):
         ct.nyquist_plot(sys, mirror_style=0.2)
 
+@pytest.mark.usefixtures("editsdefaults")
+def test_nyquist_legacy():
+    ct.use_legacy_defaults('0.9.1')
+
+    # Example that generated a warning using earlier defaults
+    s = ct.tf('s')
+    sys = (0.02 * s**3 - 0.1 * s) / (s**4 + s**3 + s**2 + 0.25 * s + 0.04)
+
+    with pytest.warns(UserWarning, match="indented contour may miss"):
+        count = ct.nyquist_plot(sys)
+
 
 if __name__ == "__main__":
     #

From 9a5759fb2809db1ca8c8f7ed42168646cc42be3a Mon Sep 17 00:00:00 2001
From: Richard Murray <murray@cds.caltech.edu>
Date: Thu, 14 Apr 2022 22:11:28 -0700
Subject: [PATCH 08/11] increase coverage; PEP8 cleanup

---
 control/freqplot.py           | 82 +++++++++++++++++++----------------
 control/tests/nyquist_test.py | 15 +++++++
 2 files changed, 60 insertions(+), 37 deletions(-)

diff --git a/control/freqplot.py b/control/freqplot.py
index 7710c0d53..006bf4830 100644
--- a/control/freqplot.py
+++ b/control/freqplot.py
@@ -149,12 +149,13 @@ def bode_plot(syslist, omega=None,
         the `deg` parameter. Default is -180 if wrap_phase is False, 0 if
         wrap_phase is True.
     wrap_phase : bool or float
-        If wrap_phase is `False`, then the phase will be unwrapped so that it
-        is continuously increasing or decreasing.  If wrap_phase is `True` the
-        phase will be restricted to the range [-180, 180) (or [:math:`-\\pi`,
-        :math:`\\pi`) radians). If `wrap_phase` is specified as a float, the
-        phase will be offset by 360 degrees if it falls below the specified
-        value. Default to `False`, set by config.defaults['freqplot.wrap_phase'].
+        If wrap_phase is `False` (default), then the phase will be unwrapped
+        so that it is continuously increasing or decreasing.  If wrap_phase is
+        `True` the phase will be restricted to the range [-180, 180) (or
+        [:math:`-\\pi`, :math:`\\pi`) radians). If `wrap_phase` is specified
+        as a float, the phase will be offset by 360 degrees if it falls below
+        the specified value. Default value is `False` and can be set using
+        config.defaults['freqplot.wrap_phase'].
 
     The default values for Bode plot configuration parameters can be reset
     using the `config.defaults` dictionary, with module name 'bode'.
@@ -573,9 +574,6 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
     return_contour : bool, optional
         If 'True', return the contour used to evaluate the Nyquist plot.
 
-    *args : :func:`matplotlib.pyplot.plot` positional properties, optional
-        Additional arguments for `matplotlib` plots (color, linestyle, etc)
-
     **kwargs : :func:`matplotlib.pyplot.plot` keyword properties, optional
         Additional keywords (passed to `matplotlib`)
 
@@ -586,15 +584,12 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
         multiple systems are given, an array of counts is returned.
 
     contour : ndarray (or list of ndarray if len(syslist) > 1)), optional
-        The contour used to create the primary Nyquist curve segment.  To
-        obtain the Nyquist curve values, evaluate system(s) along contour.
+        The contour used to create the primary Nyquist curve segment, returned
+        if `return_contour` is Tue.  To obtain the Nyquist curve values,
+        evaluate system(s) along contour.
 
     Additional Parameters
     ---------------------
-    label_freq : int, optiona
-        Label every nth frequency on the plot.  If not specified, no labels
-        are generated.
-
     arrows : int or 1D/2D array of floats, optional
         Specify the number of arrows to plot on the Nyquist curve.  If an
         integer is passed. that number of equally spaced arrows will be
@@ -630,6 +625,10 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
         imaginary axis. Portions of the Nyquist plot corresponding to indented
         portions of the contour are plotted using a different line style.
 
+    label_freq : int, optiona
+        Label every nth frequency on the plot.  If not specified, no labels
+        are generated.
+
     max_curve_magnitude : float, optional
         Restrict the maximum magnitude of the Nyquist plot to this value.
         Portions of the Nyquist plot whose magnitude is restricted are
@@ -666,6 +665,10 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
     warn_nyquist : bool, optional
         If set to 'False', turn off warnings about frequencies above Nyquist.
 
+    warn_encirclements : bool, optional
+        If set to 'False', turn off warnings about number of encirclements not
+        meeting the Nyquist criterion.
+
     Notes
     -----
     1. If a discrete time model is given, the frequency response is computed
@@ -885,22 +888,22 @@ def _parse_linestyle(style_name, allow_false=False):
                 # See if we need to indent around it
                 if abs(s - p) < indent_radius:
                     # Figure out how much to offset (simple trigonometry)
-                    offset = np.sqrt(indent_radius ** 2 - (s-p).imag ** 2) \
-                        -(s-p).real
+                    offset = np.sqrt(indent_radius ** 2 - (s - p).imag ** 2) \
+                        - (s - p).real
 
                     # Figure out which way to offset the contour point
-                    if p.real < 0 or (np.isclose(p.real, 0) \
-                            and indent_direction == 'right'):
+                    if p.real < 0 or (np.isclose(p.real, 0)
+                                      and indent_direction == 'right'):
                         # Indent to the right
                         splane_contour[i] += offset
 
-                    elif p.real > 0 or (np.isclose(p.real, 0) \
-                            and indent_direction == 'left'):
+                    elif p.real > 0 or (np.isclose(p.real, 0)
+                                        and indent_direction == 'left'):
                         # Indent to the left
                         splane_contour[i] -= offset
 
                     else:
-                        ValueError("unknown value for indent_direction")
+                        raise ValueError("unknown value for indent_direction")
 
         # change contour to z-plane if necessary
         if sys.isctime():
@@ -939,7 +942,8 @@ def _parse_linestyle(style_name, allow_false=False):
                     "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):
+            elif indent_direction == 'none' and \
+                 any(np.isclose(sys.pole().real, 0)):
                 warnings.warn(
                     "system has pure imaginary poles but indentation is"
                     " turned off; results may be meaningless",
@@ -950,14 +954,14 @@ def _parse_linestyle(style_name, allow_false=False):
 
         if plot:
             # Parse the arrows keyword
-            if isinstance(arrows, int):
+            if not arrows:
+                arrow_pos = []
+            elif isinstance(arrows, int):
                 N = arrows
                 # Space arrows out, starting midway along each "region"
                 arrow_pos = np.linspace(0.5/N, 1 + 0.5/N, N, endpoint=False)
             elif isinstance(arrows, (list, np.ndarray)):
                 arrow_pos = np.sort(np.atleast_1d(arrows))
-            elif not arrows:
-                arrow_pos = []
             else:
                 raise ValueError("unknown or unsupported arrow location")
 
@@ -1150,6 +1154,7 @@ def _add_arrows_to_line2D(
         arrows.append(p)
     return arrows
 
+
 #
 # Function to compute Nyquist curve offsets
 #
@@ -1172,13 +1177,13 @@ def _compute_curve_offset(resp, mask, max_offset):
         while i < resp.size and mask[i]:
             i += 1              # Increment the counter
             if i == resp.size:
-                break;
+                break
             # Keep track of the arclength
             arclen[i] = arclen[i-1] + np.abs(resp[i] - resp[i-1])
 
         nsegs += 0.5
         if i == resp.size:
-            break;
+            break
 
         # Save the starting offset of this segment
         seg_start = i
@@ -1187,13 +1192,13 @@ def _compute_curve_offset(resp, mask, max_offset):
         while i < resp.size and not mask[i]:
             i += 1
             if i == resp.size:  # See if we are done with this segment
-                break;
+                break
             # Keep track of the arclength
             arclen[i] = arclen[i-1] + np.abs(resp[i] - resp[i-1])
 
         nsegs += 0.5
         if i == resp.size:
-            break;
+            break
 
         # Save the ending offset of this segment
         seg_end = i
@@ -1333,7 +1338,8 @@ def singular_values_plot(syslist, omega=None,
                          *args, **kwargs):
     """Singular value plot for a system
 
-    Plots a Singular Value plot for the system over a (optional) frequency range.
+    Plots a singular value plot for the system over a (optional) frequency
+    range.
 
     Parameters
     ----------
@@ -1347,11 +1353,11 @@ def singular_values_plot(syslist, omega=None,
         Limits of the frequency vector to generate.
         If Hz=True the limits are in Hz otherwise in rad/s.
     omega_num : int
-        Number of samples to plot.
-        Default value (1000) set by config.defaults['freqplot.number_of_samples'].
+        Number of samples to plot. Default value (1000) set by
+        config.defaults['freqplot.number_of_samples'].
     dB : bool
-        If True, plot result in dB.
-        Default value (False) set by config.defaults['freqplot.dB'].
+        If True, plot result in dB. Default value (False) set by
+        config.defaults['freqplot.dB'].
     Hz : bool
         If True, plot frequency in Hz (omega must be provided in rad/sec).
         Default value (False) set by config.defaults['freqplot.Hz']
@@ -1373,7 +1379,8 @@ def singular_values_plot(syslist, omega=None,
     --------
     >>> import numpy as np
     >>> den = [75, 1]
-    >>> sys = TransferFunction([[[87.8], [-86.4]], [[108.2], [-109.6]]], [[den, den], [den, den]])
+    >>> sys = TransferFunction(
+        [[[87.8], [-86.4]], [[108.2], [-109.6]]], [[den, den], [den, den]])
     >>> omega = np.logspace(-4, 1, 1000)
     >>> sigma, omega = singular_values_plot(sys, plot=True)
     >>> singular_values_plot(sys, 0.0, plot=False)
@@ -1480,7 +1487,8 @@ def singular_values_plot(syslist, omega=None,
     # Add a grid to the plot + labeling
     if plot:
         ax_sigma.grid(grid, which='both')
-        ax_sigma.set_ylabel("Singular Values (dB)" if dB else "Singular Values")
+        ax_sigma.set_ylabel(
+            "Singular Values (dB)" if dB else "Singular Values")
         ax_sigma.set_xlabel("Frequency (Hz)" if Hz else "Frequency (rad/sec)")
 
     if len(syslist) == 1:
diff --git a/control/tests/nyquist_test.py b/control/tests/nyquist_test.py
index 8ac083535..4a9f662b3 100644
--- a/control/tests/nyquist_test.py
+++ b/control/tests/nyquist_test.py
@@ -185,6 +185,7 @@ def test_nyquist_fbs_examples():
 
 @pytest.mark.parametrize("arrows", [
     None,                       # default argument
+    False,                      # no arrows
     1, 2, 3, 4,                 # specified number of arrows
     [0.1, 0.5, 0.9],            # specify arc lengths
 ])
@@ -318,12 +319,22 @@ def test_nyquist_exceptions():
     with pytest.warns(FutureWarning, match="use `arrow_size` instead"):
         ct.nyquist_plot(sys, arrow_width=8, arrow_length=6)
 
+    # Unknown arrow keyword
+    with pytest.raises(ValueError, match="unsupported arrow location"):
+        ct.nyquist_plot(sys, arrows='uniform')
+
+    # Bad value for indent direction
+    sys = ct.tf([1], [1, 0, 1])
+    with pytest.raises(ValueError, match="unknown value for indent"):
+        ct.nyquist_plot(sys, indent_direction='up')
+
     # Discrete time system sampled above Nyquist frequency
     sys = ct.drss(2, 1, 1)
     sys.dt = 0.01
     with pytest.warns(UserWarning, match="above Nyquist"):
         ct.nyquist_plot(sys, np.logspace(-2, 3))
 
+
 def test_linestyle_checks():
     sys = ct.rss(2, 1, 1)
 
@@ -337,6 +348,10 @@ def test_linestyle_checks():
     with pytest.raises(ValueError, match="invalid 'mirror_style'"):
         ct.nyquist_plot(sys, mirror_style=0.2)
 
+    # If only one line style is given use, the default value for the other
+    # TODO: for now, just make sure the signature works; no correct check yet
+    ct.nyquist_plot(sys, primary_style=':', mirror_style='-.')
+
 @pytest.mark.usefixtures("editsdefaults")
 def test_nyquist_legacy():
     ct.use_legacy_defaults('0.9.1')

From 5dd6a70ba731983188e40e99a0519aef39b9767c Mon Sep 17 00:00:00 2001
From: Richard Murray <murray@cds.caltech.edu>
Date: Fri, 15 Apr 2022 11:08:47 -0700
Subject: [PATCH 09/11] turn off encirclement warnings for describing function
 plot

---
 control/descfcn.py            | 15 +++++++++++--
 control/freqplot.py           | 41 +++++++++++++++++++----------------
 control/tests/descfcn_test.py |  4 ++--
 3 files changed, 37 insertions(+), 23 deletions(-)

diff --git a/control/descfcn.py b/control/descfcn.py
index 2ebb18569..149db1bd2 100644
--- a/control/descfcn.py
+++ b/control/descfcn.py
@@ -199,7 +199,8 @@ def describing_function(
 
 
 def describing_function_plot(
-        H, F, A, omega=None, refine=True, label="%5.2g @ %-5.2g", **kwargs):
+        H, F, A, omega=None, refine=True, label="%5.2g @ %-5.2g",
+        warn=None, **kwargs):
     """Plot a Nyquist plot with a describing function for a nonlinear system.
 
     This function generates a Nyquist plot for a closed loop system consisting
@@ -220,6 +221,10 @@ def describing_function_plot(
     label : str, optional
         Formatting string used to label intersection points on the Nyquist
         plot.  Defaults to "%5.2g @ %-5.2g".  Set to `None` to omit labels.
+    warn : bool, optional
+        Set to True to turn on warnings generated by `nyquist_plot` or False
+        to turn off warnings.  If not set (or set to None), warnings are
+        turned off if omega is specified, otherwise they are turned on.
 
     Returns
     -------
@@ -240,9 +245,15 @@ def describing_function_plot(
     [(3.344008947853124, 1.414213099755523)]
 
     """
+    # Decide whether to turn on warnings or not
+    if warn is None:
+        # Turn warnings on unless omega was specified
+        warn = omega is None
+
     # Start by drawing a Nyquist curve
     count, contour = nyquist_plot(
-        H, omega, plot=True, return_contour=True, **kwargs)
+        H, omega, plot=True, return_contour=True,
+        warn_encirclements=warn, warn_nyquist=warn, **kwargs)
     H_omega, H_vals = contour.imag, H(contour)
 
     # Compute the describing function
diff --git a/control/freqplot.py b/control/freqplot.py
index 006bf4830..56e67e91d 100644
--- a/control/freqplot.py
+++ b/control/freqplot.py
@@ -535,9 +535,10 @@ def gen_zero_centered_series(val_min, val_max, period):
 }
 
 
-def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
-                 omega_num=None, label_freq=0, color=None,
-                 return_contour=False, warn_nyquist=True, *args, **kwargs):
+def nyquist_plot(
+        syslist, omega=None, plot=True, omega_limits=None, omega_num=None,
+        label_freq=0, color=None, return_contour=False,
+        warn_encirclements=True, warn_nyquist=True, **kwargs):
     """Nyquist plot for a system
 
     Plots a Nyquist plot for the system over a (optional) frequency range.
@@ -647,11 +648,11 @@ def nyquist_plot(syslist, omega=None, plot=True, omega_limits=None,
         determined by config.defaults['nyquist.mirror_style'].
 
     primary_style : [str, str], optional
-        Linestyles for primary image of the Nyquist curve.  The first element
-        is used for unscaled portions of the Nyquist curve, the second
-        element is used for scaled portions that are scaled (using
-        max_curve_magnitude). Default linestyle (['-', ':']) is determined by
-        config.defaults['nyquist.mirror_style'].
+        Linestyles for primary image of the Nyquist curve.  The first
+        element is used for unscaled portions of the Nyquist curve,
+        the second element is used for portions that are scaled (using
+        max_curve_magnitude). Default linestyle (['-', ':']) is
+        determined by config.defaults['nyquist.mirror_style'].
 
     start_marker : str, optional
         Matplotlib marker to use to mark the starting point of the Nyquist
@@ -839,7 +840,8 @@ def _parse_linestyle(style_name, allow_false=False):
                     p_ol = splane_poles[
                         (np.abs(splane_poles - p_cl)).argmin()]
 
-                    if abs(p_ol - p_cl) <= indent_radius:
+                    if abs(p_ol - p_cl) <= indent_radius and \
+                       warn_encirclements:
                         warnings.warn(
                             "indented contour may miss closed loop pole; "
                             "consider reducing indent_radius to be less than "
@@ -920,7 +922,8 @@ def _parse_linestyle(style_name, allow_false=False):
         count = int(np.round(encirclements, 0))
 
         # Let the user know if the count might not make sense
-        if abs(encirclements - count) > encirclement_threshold:
+        if abs(encirclements - count) > encirclement_threshold and \
+           warn_encirclements:
             warnings.warn(
                 "number of encirclements was a non-integer value; this can"
                 " happen is contour is not closed, possibly based on a"
@@ -937,13 +940,13 @@ def _parse_linestyle(style_name, allow_false=False):
             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:
+            if Z != count + P and warn_encirclements:
                 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(np.isclose(sys.pole().real, 0)):
+            elif indent_direction == 'none' and any(sys.pole().real == 0) and \
+                 warn_encirclements:
                 warnings.warn(
                     "system has pure imaginary poles but indentation is"
                     " turned off; results may be meaningless",
@@ -991,7 +994,7 @@ def _parse_linestyle(style_name, allow_false=False):
             x_reg = np.ma.masked_where(reg_mask, resp.real)
             y_reg = np.ma.masked_where(reg_mask, resp.imag)
             p = plt.plot(
-                x_reg, y_reg, primary_style[0], color=color, *args, **kwargs)
+                x_reg, y_reg, primary_style[0], color=color, **kwargs)
             c = p[0].get_color()
 
             # Figure out how much to offset the curve: the offset goes from
@@ -1005,13 +1008,13 @@ def _parse_linestyle(style_name, allow_false=False):
             y_scl = np.ma.masked_where(scale_mask, resp.imag)
             plt.plot(
                 x_scl * (1 + curve_offset), y_scl * (1 + curve_offset),
-                primary_style[1], color=c, *args, **kwargs)
+                primary_style[1], color=c, **kwargs)
 
             # Plot the primary curve (invisible) for setting arrows
             x, y = resp.real.copy(), resp.imag.copy()
             x[reg_mask] *= (1 + curve_offset[reg_mask])
             y[reg_mask] *= (1 + curve_offset[reg_mask])
-            p = plt.plot(x, y, linestyle='None', color=c, *args, **kwargs)
+            p = plt.plot(x, y, linestyle='None', color=c, **kwargs)
 
             # Add arrows
             ax = plt.gca()
@@ -1022,17 +1025,17 @@ def _parse_linestyle(style_name, allow_false=False):
             if mirror_style is not False:
                 # Plot the regular and scaled segments
                 plt.plot(
-                    x_reg, -y_reg, mirror_style[0], color=c, *args, **kwargs)
+                    x_reg, -y_reg, mirror_style[0], color=c, **kwargs)
                 plt.plot(
                     x_scl * (1 - curve_offset),
                     -y_scl * (1 - curve_offset),
-                    mirror_style[1], color=c, *args, **kwargs)
+                    mirror_style[1], color=c, **kwargs)
 
                 # Add the arrows (on top of an invisible contour)
                 x, y = resp.real.copy(), resp.imag.copy()
                 x[reg_mask] *= (1 - curve_offset[reg_mask])
                 y[reg_mask] *= (1 - curve_offset[reg_mask])
-                p = plt.plot(x, -y, linestyle='None', color=c, *args, **kwargs)
+                p = plt.plot(x, -y, linestyle='None', color=c, **kwargs)
                 _add_arrows_to_line2D(
                     ax, p[0], arrow_pos, arrowstyle=arrow_style, dir=-1)
 
diff --git a/control/tests/descfcn_test.py b/control/tests/descfcn_test.py
index 0ebfda446..796ad9034 100644
--- a/control/tests/descfcn_test.py
+++ b/control/tests/descfcn_test.py
@@ -140,7 +140,7 @@ def test_describing_function(fcn, amin, amax):
 def test_describing_function_plot():
     # Simple linear system with at most 1 intersection
     H_simple = ct.tf([1], [1, 2, 2, 1])
-    omega = np.logspace(-2, 2, 100)
+    omega = np.logspace(-1, 2, 100)
 
     # Saturation nonlinearity
     F_saturation = ct.descfcn.saturation_nonlinearity(1)
@@ -160,7 +160,7 @@ def test_describing_function_plot():
 
     # Multiple intersections
     H_multiple = H_simple * ct.tf(*ct.pade(5, 4)) * 4
-    omega = np.logspace(-2, 3, 50)
+    omega = np.logspace(-1, 3, 50)
     F_backlash = ct.descfcn.friction_backlash_nonlinearity(1)
     amp = np.linspace(0.6, 5, 50)
     xsects = ct.describing_function_plot(H_multiple, F_backlash, amp, omega)

From 5d9c7f0e890e1db7f41d5d2c4fbe93aece9cbdf9 Mon Sep 17 00:00:00 2001
From: Richard Murray <murray@cds.caltech.edu>
Date: Fri, 15 Apr 2022 17:59:07 -0700
Subject: [PATCH 10/11] rebase cleanup

---
 control/freqplot.py | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/control/freqplot.py b/control/freqplot.py
index 56e67e91d..06a85a6e1 100644
--- a/control/freqplot.py
+++ b/control/freqplot.py
@@ -937,15 +937,15 @@ def _parse_linestyle(style_name, allow_false=False):
         # 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()
+            P = (sys.poles().real > 0).sum() if indent_direction == 'right' \
+                else (sys.poles().real >= 0).sum()
+            Z = (sys.feedback().poles().real >= 0).sum()
             if Z != count + P and warn_encirclements:
                 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) and \
+            elif indent_direction == 'none' and any(sys.poles().real == 0) and \
                  warn_encirclements:
                 warnings.warn(
                     "system has pure imaginary poles but indentation is"

From 7407265496ecf03de95d31416a663afb09aa7d3c Mon Sep 17 00:00:00 2001
From: Richard Murray <murray@cds.caltech.edu>
Date: Sat, 16 Apr 2022 13:49:13 -0700
Subject: [PATCH 11/11] Fixed docstrings per roryyorke plus discrete-time
 corrections * changed return type for poles() and zeros() to complex *
 updated (missing) discrete-time tests for stability warnings * changed
 processing of near imaginary poles to pure imaginary poles

---
 control/freqplot.py           | 57 ++++++++++++++++++++++-------------
 control/statesp.py            |  8 +++--
 control/tests/nyquist_test.py | 23 +++++++++++---
 control/xferfcn.py            |  4 +--
 4 files changed, 62 insertions(+), 30 deletions(-)

diff --git a/control/freqplot.py b/control/freqplot.py
index 06a85a6e1..05ae9da55 100644
--- a/control/freqplot.py
+++ b/control/freqplot.py
@@ -644,14 +644,14 @@ def nyquist_plot(
         Linestyles for mirror image of the Nyquist curve.  The first element
         is used for unscaled portions of the Nyquist curve, the second element
         is used for portions that are scaled (using max_curve_magnitude).  If
-        `False` then omit completely.  Default linestyle (['--', '-.']) is
+        `False` then omit completely.  Default linestyle (['--', ':']) is
         determined by config.defaults['nyquist.mirror_style'].
 
     primary_style : [str, str], optional
         Linestyles for primary image of the Nyquist curve.  The first
         element is used for unscaled portions of the Nyquist curve,
         the second element is used for portions that are scaled (using
-        max_curve_magnitude). Default linestyle (['-', ':']) is
+        max_curve_magnitude).  Default linestyle (['-', '-.']) is
         determined by config.defaults['nyquist.mirror_style'].
 
     start_marker : str, optional
@@ -750,6 +750,9 @@ def _parse_linestyle(style_name, allow_false=False):
         if isinstance(style, str):
             # Only one style provided, use the default for the other
             style = [style, _nyquist_defaults['nyquist.' + style_name][1]]
+            warnings.warn(
+                "use of a single string for linestyle will be deprecated "
+                " in a future release", PendingDeprecationWarning)
         if (allow_false and style is False) or \
            (isinstance(style, list) and len(style) == 2):
             return style
@@ -765,7 +768,7 @@ def _parse_linestyle(style_name, allow_false=False):
 
     # Determine the range of frequencies to use, based on args/features
     omega, omega_range_given = _determine_omega_vector(
-        syslist, omega, omega_limits, omega_num)
+        syslist, omega, omega_limits, omega_num, feature_periphery_decades=2)
 
     # If omega was not specified explicitly, start at omega = 0
     if not omega_range_given:
@@ -790,7 +793,7 @@ def _parse_linestyle(style_name, allow_false=False):
 
         # Determine the contour used to evaluate the Nyquist curve
         if sys.isdtime(strict=True):
-            # Transform frequencies in for discrete-time systems
+            # Restrict frequencies for discrete-time systems
             nyquistfrq = math.pi / sys.dt
             if not omega_range_given:
                 # limit up to and including nyquist frequency
@@ -817,12 +820,12 @@ def _parse_linestyle(style_name, allow_false=False):
                 # because we don't need to indent for them
                 zplane_poles = sys.poles()
                 zplane_poles = zplane_poles[~np.isclose(abs(zplane_poles), 0.)]
-                splane_poles = np.log(zplane_poles)/sys.dt
+                splane_poles = np.log(zplane_poles) / sys.dt
 
                 zplane_cl_poles = sys.feedback().poles()
                 zplane_cl_poles = zplane_cl_poles[
                     ~np.isclose(abs(zplane_poles), 0.)]
-                splane_cl_poles = np.log(zplane_cl_poles)/sys.dt
+                splane_cl_poles = np.log(zplane_cl_poles) / sys.dt
 
             #
             # Check to make sure indent radius is small enough
@@ -851,8 +854,8 @@ def _parse_linestyle(style_name, allow_false=False):
             # See if we should add some frequency points near imaginary poles
             #
             for p in splane_poles:
-                # See if we need to process this pole (skip any that is on
-                # the not near or on the negative omega axis + user override)
+                # See if we need to process this pole (skip if on the negative
+                # imaginary axis or not near imaginary axis + user override)
                 if p.imag < 0 or abs(p.real) > indent_radius or \
                    omega_range_given:
                     continue
@@ -894,13 +897,13 @@ def _parse_linestyle(style_name, allow_false=False):
                         - (s - p).real
 
                     # Figure out which way to offset the contour point
-                    if p.real < 0 or (np.isclose(p.real, 0)
-                                      and indent_direction == 'right'):
+                    if p.real < 0 or (p.real == 0 and
+                                      indent_direction == 'right'):
                         # Indent to the right
                         splane_contour[i] += offset
 
-                    elif p.real > 0 or (np.isclose(p.real, 0)
-                                        and indent_direction == 'left'):
+                    elif p.real > 0 or (p.real == 0 and
+                                         indent_direction == 'left'):
                         # Indent to the left
                         splane_contour[i] -= offset
 
@@ -937,9 +940,21 @@ def _parse_linestyle(style_name, allow_false=False):
         # Nyquist criterion is actually satisfied.
         #
         if isinstance(sys, (StateSpace, TransferFunction)):
-            P = (sys.poles().real > 0).sum() if indent_direction == 'right' \
-                else (sys.poles().real >= 0).sum()
-            Z = (sys.feedback().poles().real >= 0).sum()
+            # Count the number of open/closed loop RHP poles
+            if sys.isctime():
+                if indent_direction == 'right':
+                    P = (sys.poles().real > 0).sum()
+                else:
+                    P = (sys.poles().real >= 0).sum()
+                Z = (sys.feedback().poles().real >= 0).sum()
+            else:
+                if indent_direction == 'right':
+                    P = (np.abs(sys.poles()) > 1).sum()
+                else:
+                    P = (np.abs(sys.poles()) >= 1).sum()
+                Z = (np.abs(sys.feedback().poles()) >= 1).sum()
+
+            # Check to make sure the results make sense; warn if not
             if Z != count + P and warn_encirclements:
                 warnings.warn(
                     "number of encirclements does not match Nyquist criterion;"
@@ -976,7 +991,7 @@ def _parse_linestyle(style_name, allow_false=False):
             # Find the different portions of the curve (with scaled pts marked)
             reg_mask = np.logical_or(
                 np.abs(resp) > max_curve_magnitude,
-                contour.real != 0)
+                splane_contour.real != 0)
             # reg_mask = np.logical_or(
             #     np.abs(resp.real) > max_curve_magnitude,
             #     np.abs(resp.imag) > max_curve_magnitude)
@@ -1508,7 +1523,7 @@ def singular_values_plot(syslist, omega=None,
 
 # Determine the frequency range to be used
 def _determine_omega_vector(syslist, omega_in, omega_limits, omega_num,
-                            Hz=None):
+                            Hz=None, feature_periphery_decades=None):
     """Determine the frequency range for a frequency-domain plot
     according to a standard logic.
 
@@ -1554,9 +1569,9 @@ def _determine_omega_vector(syslist, omega_in, omega_limits, omega_num,
         if omega_limits is None:
             omega_range_given = False
             # Select a default range if none is provided
-            omega_out = _default_frequency_range(syslist,
-                                                 number_of_samples=omega_num,
-                                                 Hz=Hz)
+            omega_out = _default_frequency_range(
+                syslist, number_of_samples=omega_num, Hz=Hz,
+                feature_periphery_decades=feature_periphery_decades)
         else:
             omega_limits = np.asarray(omega_limits)
             if len(omega_limits) != 2:
@@ -1640,7 +1655,7 @@ def _default_frequency_range(syslist, Hz=None, number_of_samples=None,
 
                 features_ = np.concatenate((sys.poles(), sys.zeros()))
                 # Get rid of poles and zeros on the real axis (imag==0)
-                # * origin and real < 0
+               # * origin and real < 0
                 # * at 1.: would result in omega=0. (logaritmic plot!)
                 toreplace = np.isclose(features_.imag, 0.0) & (
                                     (features_.real <= 0.) |
diff --git a/control/statesp.py b/control/statesp.py
index 58412e57a..c04174d25 100644
--- a/control/statesp.py
+++ b/control/statesp.py
@@ -987,7 +987,8 @@ def freqresp(self, omega):
     def poles(self):
         """Compute the poles of a state space system."""
 
-        return eigvals(self.A) if self.nstates else np.array([])
+        return eigvals(self.A).astype(complex) if self.nstates \
+            else np.array([])
 
     def zeros(self):
         """Compute the zeros of a state space system."""
@@ -1006,8 +1007,9 @@ def zeros(self):
             if nu == 0:
                 return np.array([])
             else:
+                # Use SciPy generalized eigenvalue fucntion
                 return sp.linalg.eigvals(out[8][0:nu, 0:nu],
-                                         out[9][0:nu, 0:nu])
+                                         out[9][0:nu, 0:nu]).astype(complex)
 
         except ImportError:  # Slycot unavailable. Fall back to scipy.
             if self.C.shape[0] != self.D.shape[1]:
@@ -1031,7 +1033,7 @@ def zeros(self):
                                            (0, self.B.shape[1])), "constant")
             return np.array([x for x in sp.linalg.eigvals(L, M,
                                                           overwrite_a=True)
-                             if not isinf(x)])
+                             if not isinf(x)], dtype=complex)
 
     # Feedback around a state space system
     def feedback(self, other=1, sign=-1):
diff --git a/control/tests/nyquist_test.py b/control/tests/nyquist_test.py
index 4a9f662b3..b1aa00577 100644
--- a/control/tests/nyquist_test.py
+++ b/control/tests/nyquist_test.py
@@ -350,7 +350,8 @@ def test_linestyle_checks():
 
     # If only one line style is given use, the default value for the other
     # TODO: for now, just make sure the signature works; no correct check yet
-    ct.nyquist_plot(sys, primary_style=':', mirror_style='-.')
+    with pytest.warns(PendingDeprecationWarning, match="single string"):
+        ct.nyquist_plot(sys, primary_style=':', mirror_style='-.')
 
 @pytest.mark.usefixtures("editsdefaults")
 def test_nyquist_legacy():
@@ -363,13 +364,17 @@ def test_nyquist_legacy():
     with pytest.warns(UserWarning, match="indented contour may miss"):
         count = ct.nyquist_plot(sys)
 
-
+def test_discrete_nyquist():
+    # Make sure we can handle discrete time systems with negative poles
+    sys = ct.tf(1, [1, -0.1], dt=1) * ct.tf(1, [1, 0.1], dt=1)
+    ct.nyquist_plot(sys)
+    
 if __name__ == "__main__":
     #
     # Interactive mode: generate plots for manual viewing
     #
-    # Running this script in python (or better ipython) will show a collection of
-    # figures that should all look OK on the screeen.
+    # Running this script in python (or better ipython) will show a
+    # collection of figures that should all look OK on the screeen.
     #
 
     # In interactive mode, turn on ipython interactive graphics
@@ -411,3 +416,13 @@ def test_nyquist_legacy():
               np.array2string(sys.poles(), precision=2, separator=','))
     count = ct.nyquist_plot(sys)
     assert _Z(sys) == count + _P(sys)
+
+    print("Discrete time systems")
+    sys = ct.c2d(sys, 0.01)
+    plt.figure()
+    plt.title("Discrete-time; poles: %s" %
+              np.array2string(sys.poles(), precision=2, separator=','))
+    count = ct.nyquist_plot(sys)
+
+    
+
diff --git a/control/xferfcn.py b/control/xferfcn.py
index 93a66ce9d..d3671c533 100644
--- a/control/xferfcn.py
+++ b/control/xferfcn.py
@@ -798,7 +798,7 @@ def poles(self):
         rts = []
         for d, o in zip(den, denorder):
             rts.extend(roots(d[:o + 1]))
-        return np.array(rts)
+        return np.array(rts).astype(complex)
 
     def zeros(self):
         """Compute the zeros of a transfer function."""
@@ -808,7 +808,7 @@ def zeros(self):
                 "for SISO systems.")
         else:
             # for now, just give zeros of a SISO tf
-            return roots(self.num[0][0])
+            return roots(self.num[0][0]).astype(complex)
 
     def feedback(self, other=1, sign=-1):
         """Feedback interconnection between two LTI objects."""