Skip to content

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.

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
merged 1 commit into from
Apr 27, 2022

Conversation

roryyorke
Copy link
Contributor

linfnorm requires Slycot routine ab13dd, which does the actual work.

Added tests for SISO and MIMO continuous- and discrete-time. MIMO
reference values are generated by linfnorm itself, but results have
been reviewed graphically for correctness.

The MIMO tests rely, unfortunately, on circular logic: the very function being tested has been used to generate test data. I could scour the various textbooks I have for independent examples, but I'm more-or-less assuming that SLICOT's AB13DD is correct, and the tests are to make sure no gross errors have occurred in wrapping. I've tried to cover the various possible cases of gpeak finite (but non-zero) and infinite, and fpeak 0, finite and non-zero, and infinite. I just realized I haven't checked for gpeak = 0; I'll add that (system would have to be non-minimal, or be a static gain of 0; I'll try the former first).

Below are the MIMO test cases for review; everything looks OK to me. The left-hand figures are continuous-time, the right-hand discrete time.

While doing this rather hacky test-generation, I found that StateSpace's repr is not eval-able; this is due, in turn, to numpy.ndarray's repr not being eval-able. I don't know why that is, but I monkey-patched StateSpace.__repr__ to address this. This new repr is eval-able, and I think looks a little better too (look in linfnorm_mimo_testcases.py for examples).

gdamped
ghigh
gint
gstatic
gun
gunder

@roryyorke
Copy link
Contributor Author

I don't know how to handle those np.matrix "PendingDeprecation" warnings-as-errors. I see matrixfilter and matrixerrorfilter marks in conftest.py, but I'm not sure when they're supposed to be used.

@roryyorke
Copy link
Contributor Author

I am still a bit lost with the PendingDeprecation error from Numpy, which I can't reproduce locally. I developed this PR in a conda environment with Python 3.9.2, pytest 6.2.3 (which I see is a bit old), numpy 1.20.2, scipy 1.6.2, and Slycot 0.4.0.0. python-control tests pass in this environment.

I just created a brand-new conda environment with Python 3.7.12, Pytest 7.1.2, Numpy 1.21.6, Scipy 1.7.3, and Slycot 0.4.0.0. python-control tests all pass in this environment.

I'm running on a 64-bit Intel i5 Ubuntu 20.04 system.

@bnavigator
Copy link
Contributor

Your environment should be perfectly fine. I guess you can reproduce locally when you set export PYTHON_CONTROL_ARRAY_AND_MATRIX=1.

The critical part is not to create StateSpace or TransferFunction objects outside of the actual test functions. They are currently created at import time of the module, which is during the test collection.The filter fixtures don't apply yet. You have to move the systems into reference fixtures.

@roryyorke
Copy link
Contributor Author

You have to move the systems into reference fixtures.

Thanks, I'll investigate that.

Not being able to reproduce the error locally makes this tricky.

@coveralls
Copy link

coveralls commented Apr 24, 2022

Coverage Status

Coverage increased (+0.006%) to 94.312% when pulling 02b3451 on roryyorke:rory/linfnorm into 46b3c53 on python-control:master.

@roryyorke
Copy link
Contributor Author

export PYTHON_CONTROL_ARRAY_AND_MATRIX=1

Oops, missed that. At least it's passing now.

Copy link
Member

@murrayrm murrayrm left a comment

Choose a reason for hiding this comment

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

Some small comments.

The method of generating test cases is a bit nonstandard. Do we need to check correctness of the algorithm itself (versus the processing of arguments, etc)? Seems like we can assume slycot is working and just do a few simple unit tests??

-------
gpeak : non-negative scalar
L-infinity norm
fpeak : non-negatie scalar
Copy link
Member

Choose a reason for hiding this comment

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

negatie → negative

Comment on lines 1921 to 1922
For systems with an undefined time base (i.e., sys.dt is None), a
ValueError is raised.
Copy link
Member

Choose a reason for hiding this comment

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

Normally we allow dt=None to be treat as either continuous or discrete time, with default being continuous. Should we do that here? This is not a major issue since dt=0 by default, but there are situations where a constant system can be converted to a state space or transfer function with dt=None (since it works as either discrete or continuous).

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I'll drop the ValueError, and use .isctime() as suggested below; this will treat these as continuous-time.

One case where dt=None by default is static gains (this is from observation, not from reading the code); and, in this case, it doesn't matter what the time base is.

Comment on lines 1943 to 1945
if sys.dt is None:
raise ValueError('unspecified time base')
dico = 'C' if sys.dt==0 else 'D'
Copy link
Member

Choose a reason for hiding this comment

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

Alternatively:

dico = 'C' if sys.isctime() else 'D'

Comment on lines +1935 to +1942
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
Copy link
Member

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.

@roryyorke
Copy link
Contributor Author

The method of generating test cases is a bit nonstandard. Do we need to check correctness of the algorithm itself (versus the processing of arguments, etc)? Seems like we can assume slycot is working and just do a few simple unit tests??

Fair enough, I did get a bit carried away.

We do need to test at least for:

  • static gains (ab13dd doesn't handle these)
  • calling ab13dd with correct DICO value for continuous or discrete-time
  • scaling of frequency for discrete-time systems

What do you think about keeping one continuous-time SISO test (a slightly underdamped system?), an equivalent SISO discrete-time test, and one continuous-time MIMO test (probably duplicating the SISO test as a 2x2 diagonal system) ? Also one static-gain test.

@roryyorke
Copy link
Contributor Author

All points addressed. If the tests pass, I'll squash the PR.

linfnorm requires Slycot routine ab13dd, which does the actual work.

Added tests to check that wrapping of ab13dd is correct, that static
systems are handled, and that discrete-time system frequencies are
scaled correctly.
Copy link
Contributor

@bnavigator bnavigator left a comment

Choose a reason for hiding this comment

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

LGTM

@murrayrm murrayrm merged commit dd95f35 into python-control:master Apr 27, 2022
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

4 participants