GV

GeezerVison

The World Wide Web, restored to 800 pixels of glory.

Home|About|Gaming cabinet|Archives|How it works|Starter sites

On this site

  • Sign in
  • ActionsCI/CD and automation
  • CodespacesCloud development environments
  • All Features
  • Documentation
  • Training
  • Certifications
  • Build with Copilot
  • Customize Copilot
  • Integrations
  • Trending repos
  • GitHub Sponsors
  • Developer articles
  • Enterprises

This view

  • Original page
  • New address

GeezerVison is showing github.com — ENH: add `linform` to compute linear system L-infinity norm by roryyorke · Pull Request #729 · python-control/python-control

ENH: add linform to compute linear system L-infinity norm #729

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.

Sign up for GitHub

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

Already on GitHub? Sign in to your account

Merged
murrayrm merged 1 commit into python-control:master from roryyorke:rory/linfnorm
Apr 27, 2022

Uh oh!

There was an error while loading. Please reload this page.

Merged

ENH: add linform to compute linear system L-infinity norm #729

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
02b3451
ENH: add `linform` to compute linear system L-infinity norm
roryyorke Apr 24, 2022
File filter
Filter by extension

Clear filters
Conversations
Failed to load comments.
Loading

Uh oh!

There was an error while loading. Please reload this page.

Jump to
Jump to file
Failed to load files.
Loading

Uh oh!

There was an error while loading. Please reload this page.

Diff view
Diff view

There are no files selected for viewing

66 changes: 65 additions & 1 deletion 66 control/statesp.py
Show comments
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -55,6 +55,7 @@
from numpy.linalg import solve, eigvals, matrix_rank
from numpy.linalg.linalg import LinAlgError
import scipy as sp
import scipy.linalg
from scipy.signal import cont2discrete
from scipy.signal import StateSpace as signalStateSpace
from warnings import warn
Expand All @@ -63,7 +64,12 @@
from . import config
from copy import deepcopy

__all__ = ['StateSpace', 'tf2ss', 'ssdata']
try:
from slycot import ab13dd
except ImportError:
ab13dd = None

__all__ = ['StateSpace', 'tf2ss', 'ssdata', 'linfnorm']

# Define module default parameter values
_statesp_defaults = {
Expand Down Expand Up @@ -1888,3 +1894,61 @@ def ssdata(sys):
"""
ss = _convert_to_statespace(sys)
return ss.A, ss.B, ss.C, ss.D


def linfnorm(sys, tol=1e-10):
"""L-infinity norm of a linear system

Parameters
----------
sys : LTI (StateSpace or TransferFunction)
system to evalute L-infinity norm of
tol : real scalar
tolerance on norm estimate

Returns
-------
gpeak : non-negative scalar
L-infinity norm
fpeak : non-negative scalar
Frequency, in rad/s, at which gpeak occurs

For stable systems, the L-infinity and H-infinity norms are equal;
for unstable systems, the H-infinity norm is infinite, while the
L-infinity norm is finite if the system has no poles on the
imaginary axis.

See also
--------
slycot.ab13dd : the Slycot routine linfnorm that does the calculation
"""

if ab13dd is None:
raise ControlSlycot("Can't find slycot module 'ab13dd'")

a, b, c, d = ssdata(_convert_to_statespace(sys))
e = np.eye(a.shape[0])

n = a.shape[0]
m = b.shape[1]
p = c.shape[0]

if n == 0:
# ab13dd doesn't accept empty A, B, C, D;
# static gain case is easy enough to compute
gpeak = scipy.linalg.svdvals(d)[0]
# max svd is constant with freq; arbitrarily choose 0 as peak
fpeak = 0
return gpeak, fpeak
Comment thread
Comment on lines +1936 to +1942

@murrayrm murrayrm Apr 25, 2022

Copy link
Copy Markdown
Member

There was a problem hiding this comment.

Choose a reason for hiding this comment

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

Note that this will return a value even if sys.dt == None, contrary to the documentation.

Sorry, something went wrong.

Uh oh!

There was an error while loading. Please reload this page.

All reactions

dico = 'C' if sys.isctime() else 'D'
jobe = 'I'
equil = 'S'
jobd = 'Z' if all(0 == d.flat) else 'D'

gpeak, fpeak = ab13dd(dico, jobe, equil, jobd, n, m, p, a, e, b, c, d, tol)

if dico=='D':
fpeak /= sys.dt

return gpeak, fpeak
53 changes: 52 additions & 1 deletion 53 control/tests/statesp_test.py
Show comments
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -19,13 +19,15 @@
from control.dtime import sample_system
from control.lti import evalfr
from control.statesp import StateSpace, _convert_to_statespace, tf2ss, \
_statesp_defaults, _rss_generate
_statesp_defaults, _rss_generate, linfnorm
from control.iosys import ss, rss, drss
from control.tests.conftest import ismatarrayout, slycotonly
from control.xferfcn import TransferFunction, ss2tf


from .conftest import editsdefaults


class TestStateSpace:
"""Tests for the StateSpace class."""

Expand Down Expand Up @@ -1097,3 +1099,52 @@ def test_latex_repr_testsize(editsdefaults):

gstatic = ss([], [], [], 1)
assert gstatic._repr_latex_() is None


class TestLinfnorm:
# these are simple tests; we assume ab13dd is correct
# python-control specific behaviour is:
# - checking for continuous- and discrete-time
# - scaling fpeak for discrete-time
# - handling static gains

# the underdamped gpeak and fpeak are found from
# gpeak = 1/(2*zeta*(1-zeta**2)**0.5)
# fpeak = wn*(1-2*zeta**2)**0.5
@pytest.fixture(params=[
('static', ct.tf, ([1.23],[1]), 1.23, 0),
('underdamped', ct.tf, ([100],[1, 2*0.5*10, 100]), 1.1547005, 7.0710678),
])
def ct_siso(self, request):
name, systype, sysargs, refgpeak, reffpeak = request.param
return systype(*sysargs), refgpeak, reffpeak

@pytest.fixture(params=[
('underdamped', ct.tf, ([100],[1, 2*0.5*10, 100]), 1e-4, 1.1547005, 7.0710678),
])
def dt_siso(self, request):
name, systype, sysargs, dt, refgpeak, reffpeak = request.param
return ct.c2d(systype(*sysargs), dt), refgpeak, reffpeak

@slycotonly
def test_linfnorm_ct_siso(self, ct_siso):
sys, refgpeak, reffpeak = ct_siso
gpeak, fpeak = linfnorm(sys)
np.testing.assert_allclose(gpeak, refgpeak)
np.testing.assert_allclose(fpeak, reffpeak)

@slycotonly
def test_linfnorm_dt_siso(self, dt_siso):
sys, refgpeak, reffpeak = dt_siso
gpeak, fpeak = linfnorm(sys)
# c2d pole-mapping has round-off
np.testing.assert_allclose(gpeak, refgpeak)
np.testing.assert_allclose(fpeak, reffpeak)

@slycotonly
def test_linfnorm_ct_mimo(self, ct_siso):
siso, refgpeak, reffpeak = ct_siso
sys = ct.append(siso, siso)
gpeak, fpeak = linfnorm(sys)
np.testing.assert_allclose(gpeak, refgpeak)
np.testing.assert_allclose(fpeak, reffpeak)
Add this suggestion to a batch that can be applied as a single commit. This suggestion is invalid because no changes were made to the code. Suggestions cannot be applied while the pull request is closed. Suggestions cannot be applied while viewing a subset of changes. Only one suggestion per line can be applied in a batch. Add this suggestion to a batch that can be applied as a single commit. Applying suggestions on deleted lines is not supported. You must change the existing code in this line in order to create a valid suggestion. Outdated suggestions cannot be applied. This suggestion has been applied or marked resolved. Suggestions cannot be applied from pending reviews. Suggestions cannot be applied on multi-line comments. Suggestions cannot be applied while the pull request is queued to merge. Suggestion cannot be applied right now. Please check back later.
You are visitor 01402503 — GeezerVison does not store the pages you look up.
Best viewed with Netscape Navigator 4.7 or Internet Explorer 6 at 800×600. No cookies required.