Movatterモバイル変換

[0]ホーム
URL:

Skip to content

Navigation Menu

Sign in
Appearance settings

Search code, repositories, users, issues, pull requests...

Provide feedback

We read every piece of feedback, and take your input very seriously.

Saved searches

Use saved searches to filter your results more quickly

 Sign up
Appearance settings

Reimplementation of 2D phase plots#980

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to ourterms of service andprivacy statement. We’ll occasionally send you account related emails.

Already on GitHub?Sign in to your account

Merged
murrayrm merged 4 commits intopython-control:mainfrommurrayrm:phaseplots_01Jan2024
Mar 31, 2024

Conversation

murrayrm
Copy link
Member

This PR reimplements the phase plot functionality of python-control by creating a new functionphase_plane_plot and creating a modulephaseplot that provides more specialized functionality. The legacyphase_plot function is still available (with a deprecation warning).

From the documentation:

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

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

phaseplot-dampedosc-default

By default, the plot includes streamlines generated from starting points on limits of the plot, with arrows showing the flow of the system, as well as any equilibrium points for the system. A variety of options are available to modify the information that is plotted, including plotting a grid of vectors instead of streamlines 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*pi, 2*pi, -2, 2], 5,    gridtype='meshgrid', gridspec=[5, 8], arrows=3,    plot_equilpoints={'gridspec': [12, 9]},    params={'m': 1, 'l': 1, 'b': 0.2, 'g': 1})plt.xlabel(r"$\theta$ [rad]")plt.ylabel(r"$\dot\theta$ [rad/sec]")

phaseplot-invpend-meshgrid

This figure shows several features of more complex phase plane plots: multiple equilibrium points are shown, with saddle points showing separatrices, and streamlines generated along a 5x8 mesh of initial conditions. At each mesh point, a streamline is created that goes 5 time units forward and backward in time. A separate grid specification is used to find equilibrium points and separatrices (since the course grid spacing of 5x8 does not find all possible equilibrium points). Together, the multiple features in the phase plane plot give a good global picture of the topological strucrure of solutions of the dynamical system.

Phase plots can be buit up by hand using a variety of helper functions that are part of thecontrol.phaseplot (pp) module:

import control.phaseplot as ppdef 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)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*pi, arrows=6, color='b')plt.gca().set_aspect('equal')

phaseplot-oscillator-helpers

The following helper functions are available:

  • pp.equilpoints
  • pp.separatrices
  • pp.streamlines
  • pp.vectorfield

Thephase_plane_plot function calls these helper functions based on the options it is passed.

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

@coveralls
Copy link

coveralls commentedMar 28, 2024•
edited
Loading

Coverage Status

coverage: 94.467% (+0.05%) from 94.422%
when pulling748e36a on murrayrm:phaseplots_01Jan2024
into7a70be1 on python-control:main.

@murrayrmmurrayrm added this to the0.10.0 milestoneMar 28, 2024
@@ -35,7 +35,7 @@
# Author: Richard M. Murray
# Date: 1 Jul 2019

r"""The :mod:`control.flatsys`package contains a set of classes and functions
r"""The :mod:`control.flatsys`module contains a set of classes and functions
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

The original is correct:flatsys is a package according to standard terms:https://docs.python.org/3/glossary.html#term-package

import numpy as np
import matplotlib.pyplot as mpl
"""The :mod:`control.phaseplot` module contains functions for generating 2D
phase plots.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Convention is for summary to be 1 line (https://peps.python.org/pep-0257/#multi-line-docstrings).

@@ -0,0 +1,215 @@
# phase_portraits.py - phase portrait examples
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Suggested change
#phase_portraits.py - phase portrait examples
#phase_plane_plots.py - phase portrait examples

@@ -29,6 +30,9 @@ freqplot-siso_bode-default.png: ../control/tests/freqplot_test.py
rlocus-siso_ctime-default.png: ../control/tests/rlocus_test.py
PYTHONPATH=.. python $<

phaseplot-dampedosc-default.png: ../control/tests/phaseplot_test.py
PYTHONPATH=.. python $<
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Should the other images generated by phaseplot_test.py (phaseplot-oscillator-helpers.png,phaseplot-invpend-meshgrid.png) be included here? It does not seem important because the PNG files are being committed, but I note this here in case you wanted to.

Copy link
MemberAuthor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Agree that they could be listed, but I didn't list them all out since triggering on the one file will update all files => if the test generation script is updated, all PNG files are updated.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Good point! I will complete review tomorrow.

Copy link
Contributor

@sawyerbfullersawyerbfuller left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Fixes some very outdated code. Very nice!

return new_kwargs

# Create array for storing outputs
out = np.empty(3, dtype=object)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Docstring says this is a list. Seems like a list might be better/simpler here than a numpy array, e.g.out = [[], None, None]


Parameters
----------
sys : NonlinearIOSystem or callable(x, t, ...)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

What is motivation to maintain the old convention of (x, t, ...) argument order, rather than (t, x, ...) when called with a callable? Backwards compatibility? I am not aware of anywhere else we or scipy adheres to that convention anymore, but I might be missing something.

vfdata = np.zeros((points.shape[0], 4))
for i, x in enumerate(points):
vfdata[i, :2] = x
vfdata[i, 2:] = sys.dynamics(0, x, 0, params=params)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

I think this wouldn't work ifsys is just a function - is the callable intended to be a function if it is not a system?

Parameters
----------
sys : NonlinearIOSystem or callable(x, t, ...)
Function used to generate phase plane data.
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Relatedly, should there be a statement here that only two-state systems can be accepted? (and a test that gives a helpful error?)

@murrayrm
Copy link
MemberAuthor

In addition to adressing comments, functionality for handling callables (versus systems) was missing; now added.

# Create reverse time system, if needed
if dir != 'forward':
revsys = NonlinearIOSystem(
lambda t, x, u, params: -np.array(sys.updfcn(t, x, u, params)),
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Suggested change
lambdat,x,u,params:-np.array(sys.updfcn(t,x,u,params)),
lambdat,x,u,params:-np.asarray(sys.updfcn(t,x,u,params)),

This won't re-create the array if it is already one.

relatedly, and this is more of a question: I am confused about the use of each of these methods, all of which do roughly the same thing. Is this correct:_rhs fast and low-level,updfcn just a reference to the function that was passed in, anddynamics: user-friendly function that may not be as fast? Just checking that this is the right one (speed is kind of important in this application)

Copy link
MemberAuthor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Correct on the various versions. In general:

  • _updfcn() should not be called directly (though OK here, because we are creating a new system).
  • _rhs() should be called inside loops (with a call to _update_params() before the loop starts, which sets the parameter values)
  • dynamics() is the user-callable function (which calls _update_params and then _rhs)

Copy link
Contributor

@sawyerbfullersawyerbfuller left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

LGTM

Copy link
Member

@slivingstonslivingston left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Great work! I found a few more misprints, but nothing critical. Locally I merged this into the current tip ofmain branch and successfully ran the tests and example code, and built the documentation.

# RMM, 25 Mar 2024
#
# This file contains a number of examples of phase plane plots generated
# using the phaseplot module. Most of these figures lines up with examples
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Suggested change
# using the phaseplot module. Most of these figureslines up with examples
# using the phaseplot module. Most of these figuresline up with examples

phase plots. The base function for creating phase plane portraits is
:func:`~control.phase_plane_plot`, which generates a phase plane portrait
for a 2 state I/O system (with no inputs). In addition, several other
functions are available to creat ecustomize phase plane plots:
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Suggested change
functionsareavailabletocreatecustomizephaseplaneplots:
functionsareavailabletocustomizephaseplaneplots:

Parameters
----------
sys : NonlinearIOSystem or callable(t, x, ...)
I/O systems or function used to generate phase plane data. If a
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Suggested change
I/Osystemsorfunctionusedtogeneratephaseplanedata.Ifa
I/Osystemorfunctionusedtogeneratephaseplanedata.Ifa

The :func:`~control.phase_plane_plot` function calls these helper functions
based on the options it is passed.

Note that unlike other plotting function, phase plane plots do not involve
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Suggested change
Note that unlike other plottingfunction, phase plane plots do not involve
Note that unlike other plottingfunctions, phase plane plots do not involve

plot_separatrices : bool or dict
If `True` (default) then plot separatrices starting from each
equilibrium point. If set to a dict, pass on the key-value pairs
in the dict as keywords to :func:`~control.phaseplot.streamlines`.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Suggested change
inthedictaskeywordsto :func:`~control.phaseplot.streamlines`.
inthedictaskeywordsto :func:`~control.phaseplot.separatrices`.


import numpy as np
import matplotlib.pyplot as mpl
# Create a system form a callable
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Suggested change
# Create a systemform a callable
# Create a systemfrom a callable

Comment on lines 792 to 793
# if not isinstance(pointdata, (list, tuple)) or len(pointdata) != 4:
# raise ValueError("invalid grid data specification")
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others.Learn more.

Suggested change
# if not isinstance(pointdata, (list, tuple)) or len(pointdata) != 4:
# raise ValueError("invalid grid data specification")

This can be deleted, or the check in the comment can be expanded to handle case ofisinstance(pointdata, np.ndarray).

@murrayrm
Copy link
MemberAuthor

Thanks for the careful reviews,@sawyerbfuller and@slivingston! Will merge as soon as checks finish up (just in case).

@murrayrmmurrayrm merged commit8983e39 intopython-control:mainMar 31, 2024
@murrayrmmurrayrm deleted the phaseplots_01Jan2024 branchMarch 31, 2024 05:09
Sign up for freeto join this conversation on GitHub. Already have an account?Sign in to comment
Reviewers

@slivingstonslivingstonslivingston approved these changes

@sawyerbfullersawyerbfullersawyerbfuller approved these changes

Assignees
No one assigned
Labels
None yet
Projects
None yet
Milestone
0.10.0
4 participants
@murrayrm@coveralls@slivingston@sawyerbfuller
[8]ページ先頭
©2009-2025 Movatter.jp