Skip to content
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

Dgebal application order #122

Merged
merged 5 commits into from
May 16, 2020

Conversation

repagh
Copy link
Member

@repagh repagh commented May 12, 2020

Should address issue #11

Copy link
Collaborator

@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.

I agree with the fix in the Fortran code. Seems consistent with the DGEBAL documentation now. Good thinking to look for other possible uses of DGEBAL.

@@ -11,6 +11,7 @@
from numpy.testing import assert_almost_equal



Copy link
Collaborator

Choose a reason for hiding this comment

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

E303, too many blank lines

Suggested change

Copy link
Member Author

Choose a reason for hiding this comment

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

Should be fixed in the commit.

Comment on lines 52 to 66
def test_tb05ad_ag(self):
"""
Test that tb05ad with job 'AG' computes the correct
frequency response.
"""
for key in CASES:
sys = CASES[key]
self.check_tb05ad_AG_NG(sys, 10*1j, 'AG')

def test_tb05ad_ag_fixed_bug_no11(self):
""" Test tb05ad and job 'AG' (i.e., balancing enabled).

Failed on certain A matrices before, issue #11.
"""
self.check_tb05ad_AG_NG(CASES['fail1'], 10*1j, 'AG')
Copy link
Collaborator

Choose a reason for hiding this comment

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

This just runs CASES['fail1'] twice. Probably the name 'fail1' is not appropriate anymore.

Suggested change
def test_tb05ad_ag(self):
"""
Test that tb05ad with job 'AG' computes the correct
frequency response.
"""
for key in CASES:
sys = CASES[key]
self.check_tb05ad_AG_NG(sys, 10*1j, 'AG')
def test_tb05ad_ag_fixed_bug_no11(self):
""" Test tb05ad and job 'AG' (i.e., balancing enabled).
Failed on certain A matrices before, issue #11.
"""
self.check_tb05ad_AG_NG(CASES['fail1'], 10*1j, 'AG')
def test_tb05ad_ag(self):
"""
Test that tb05ad with job 'AG' computes the correct
frequency response with balancing A first.
Failed on certain A matrices before, issue #11.
"""
for key in CASES:
sys = CASES[key]
self.check_tb05ad_AG_NG(sys, 10*1j, 'AG')

Copy link
Member Author

Choose a reason for hiding this comment

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

renamed, and removed the repeat

def test_tb05ad_balance(self):
"""Specifically check for the balancing output.

Based on https://rdrr.io/rforge/expm/man/balance.html
Copy link
Collaborator

Choose a reason for hiding this comment

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

According to that link, this case has SCALE=[1 1 1 3], ILO=2, IHI=3 and thus only swaps 4 with 3 and 1 with 1. Does not really check the correct order of transformations done, does it?

To check the newly corrected transformation in TB05AD, Bt, Ct need to be checked for correct order as well.

Copy link
Member Author

Choose a reason for hiding this comment

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

I created a new test, based on a random matrix, tweaked until that one showed permutation

Copy link
Collaborator

Choose a reason for hiding this comment

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

Could this be pinned to a fixed found matrix? Having random data in tests that could fail at times and does not at other times will make debugging hard.

Copy link
Member Author

Choose a reason for hiding this comment

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

we're still covered by the random seed at the start of the file

Copy link
Collaborator

Choose a reason for hiding this comment

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

Ah okay, I didn't see that.

Comment on lines 234 to 237
# and verify that eigenvalues of both A matrices are close
assert_almost_equal(np.dot(np.dot(Ct, At), Bt),
np.dot(np.dot(C, A), B))
assert_almost_equal(eig(At)[0], eig(A)[0])
Copy link
Collaborator

Choose a reason for hiding this comment

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

Is the order and scaling of eigenvalues guaranteed? Could run into assertion failures on some build systems. Or do scipy.linalg.eig() and MB05AD call the same LAPACK code and thus always yield the same results?

Copy link
Member Author

Choose a reason for hiding this comment

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

I'll put in some sorting.

Copy link
Collaborator

Choose a reason for hiding this comment

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

What about scaling, including different signs? I had to deal with that in test_mb05md() (#79)

Copy link
Member Author

Choose a reason for hiding this comment

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

Now you're confusing me. Do you mean scaling of the eigenvalues? But then the frequency should also be scaled to calculate the response, and I don't see anything on that in the TB05AD.f code. Eigenvectors may be differently scaled, but I did not do a test on that.

Or do you mean the scaling in the matrix_balance / dgebal step. I now added some code for that, got to the same testing case as before, on the 10th matrix tried.

Copy link
Collaborator

Choose a reason for hiding this comment

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

Okay strike that. Brainfart by me. Mixed up eigenvalue properties with eigenvectors.

@repagh repagh merged commit b23ceb8 into python-control:master May 16, 2020
@repagh repagh mentioned this pull request May 21, 2020
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.

2 participants