Phase Plane Plots

Insight into nonlinear systems can often be obtained by looking at phase plane diagrams. The phase_plane_plot() function allows the creation of a 2-dimensional phase plane diagram for a system. This functionality is supported by a set of mapping functions that are part of the phaseplot module.

The default method for generating a phase plane plot is to provide a 2D dynamical system along with a range of coordinates in phase space:

def sys_update(t, x, u, params):
    return np.array([[0, 1], [-1, -1]]) @ x
sys = ct.nlsys(
    sys_update, states=['position', 'velocity'],
    inputs=0, name='damped oscillator')
axis_limits = [-1, 1, -1, 1]
ct.phase_plane_plot(sys, axis_limits)

By default the plot includes streamlines infered from function values on a grid, equilibrium points and separatrices if they exist. A variety of options are available to modify the information that is plotted, including plotting a grid of vectors instead of streamlines, plotting streamlines from arbitrary starting points and turning on and off various features of the plot.

To illustrate some of these possibilities, consider a phase plane plot for an inverted pendulum system, which is created using a mesh grid:

def invpend_update(t, x, u, params):
    m, l, b, g = params['m'], params['l'], params['b'], params['g']
    return [x[1], -b/m * x[1] + (g * l / m) * np.sin(x[0]) + u[0]/m]
invpend = ct.nlsys(invpend_update, states=2, inputs=1, name='invpend')

ct.phase_plane_plot(
    invpend, [-2 * np.pi, 2 * np.pi, -2, 2],
    params={'m': 1, 'l': 1, 'b': 0.2, 'g': 1})
plt.xlabel(r"$\theta$ [rad]")
plt.ylabel(r"$\dot\theta$ [rad/sec]")

This figure shows several features of more complex phase plane plots: multiple equilibrium points are shown, with saddle points showing separatrices, and streamlines generated generated from a rectangular 25x25 grid (default) of function evaluations. Together, the multiple features in the phase plane plot give a good global picture of the topological structure of solutions of the dynamical system.

Phase plots can be built up by hand using a variety of helper functions that are part of the phaseplot (pp) module. For more precise control, the streamlines can also generated by integrating the system forwards or backwards in time from a set of initial conditions. The initial conditions can be chosen on a rectangular grid, rectangual boundary, circle or from an arbitrary set of points.

import control.phaseplot as pp

def oscillator_update(t, x, u, params):
    return [x[1] + x[0] * (1 - x[0]**2 - x[1]**2),
            -x[0] + x[1] * (1 - x[0]**2 - x[1]**2)]
oscillator = ct.nlsys(
    oscillator_update, states=2, inputs=0, name='nonlinear oscillator')

ct.phase_plane_plot(oscillator, [-1.5, 1.5, -1.5, 1.5], 0.9,
    plot_streamlines=True)
pp.streamlines(
    oscillator, np.array([[0, 0]]), 1.5,
    gridtype='circlegrid', gridspec=[0.5, 6], dir='both')
pp.streamlines(
    oscillator, np.array([[1, 0]]), 2 * np.pi, arrows=6, color='b')
plt.gca().set_aspect('equal')

The following helper functions are available:

phaseplot.equilpoints(sys, pointdata[, ...])	Plot equilibrium points in the phase plane.
phaseplot.separatrices(sys, pointdata[, ...])	Plot separatrices in the phase plane.
phaseplot.streamlines(sys, pointdata[, ...])	Plot stream lines in the phase plane.
phaseplot.streamplot(sys, pointdata[, ...])	Plot streamlines in the phase plane.
phaseplot.vectorfield(sys, pointdata[, ...])	Plot a vector field in the phase plane.

The phase_plane_plot() function calls these helper functions based on the options it is passed.

Note that unlike other plotting functions, phase plane plots do not involve computing a response and then plotting the result via a plot() method. Instead, the plot is generated directly be a call to the phase_plane_plot() function (or one of the phaseplot helper functions).