[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: modified residues() for matlab compatibility
From: |
Ben Abbott |
Subject: |
Re: modified residues() for matlab compatibility |
Date: |
Mon, 8 Oct 2007 22:32:56 -0400 |
On Oct 8, 2007, at 3:41 PM, John W. Eaton wrote:
On 6-Oct-2007, Ben Abbott wrote:
| I've attached a modified residue.m that includes modifications to
the
| doc-string and the addition of test scripts. In addition, the
| multiplicity vector, "e", which the prior change dropped has been
| reinstated.
OK, I applied this patch, but had to do some of it by hand since there
had been intervening formatting changes. Would you please check the
current CVS version and make sure that it is correct? In particular,
I am seeing the failure:
octave:1> test residue
***** test
b = [1, 0, 1];
a = [1, 0, 18, 0, 81];
[r, p, k, e] = residue(b, a);
assert ((abs (54*r - [-5i; 12; 5i; 12]) < 1e-6
&& abs (p - [3i; 3i; -3i; -3i]) < 1e-7
&& isempty (k)
&& e == [1; 2; 1; 2]));
[br, ar] = residue (r, p, k);
assert ((abs (br - b) < 1e-12
&& abs (ar - a) < 1e-12));
!!!!! test failed
error: assert ((abs (54 * r - [-5i; 12; 5i; 12]) < 1e-6 && abs (p
- [3i; 3i; -3i; -3i]) < 1e-7 && isempty (k) && e == [1; 2; 1; 2]))
failed
octave:2> b = [1, 0, 1];
octave:3> a = [1, 0, 18, 0, 81];
octave:4> [r, p, k, e] = residue(b, a);
octave:5> abs (54*r - [-5i; 12; 5i; 12])
ans =
1.0000e+01
4.4460e-08
1.0000e+01
4.4460e-08
and I'm not sure whether it is the test that is bad, or residue that
is computing the wrong result.
| The test scripts include the examples give in the doc-string as well
| as the original error which recently prompted much discussion and
| lead to my involvement. The accuracy of the residues may have been
| impacted. I had to lower the increase the permitted error in the
test
| scripts in order for the scripts to pass ... at least relative to
the
| test script in your prior email. Specifically
|
| -%! assert ((abs (r - [-2; 7; 3]) < 1e-7
| +%! assert ((abs (r - [-2; 7; 3]) < 1e-5
|
| It occurs to me that I need to get set up with CVS so I can see the
| present state of things. That way I can send patches instead of
| sources ... unfortunately the web-site, http://www.octave.org/
| download.html, is down. I'll get to that task when it's back up.
I should be up again.
| In any event, I've attached a new residue.m as well as the
equivalent
| patch (which may have a unix/dos EOL problem). I've also attached a
| ChangeLog.
|
| Some questions,
|
| (1) Is there a naming convention for patches?
No the name doesn't matter, but it is helpful if you send them as
plain text attachments so I can apply them directly from Emacs instead
of having to save them to a file.
| (2) Which variety of EOL is desired, or does it matter?
The files should have Unix line endings (LF only), but I don't think
that matters if you send text/plain attachments.
thanks,
jwe
John,
I did not get the same error. However, I did get an error from the
final test script placed on the bottom of the file. Since I was
unable to connect to the CVS over the weekend and I did now know if
there is a protocol for test-scripts, I'd placed my testing between
the doc-strip and remainder of the function.
Since having access to the CVS, I've moved the scripts to the bottom
and reconciled the duplicate.
However, since this is my 1st time using the "test" functionality, it
might be wise to double check my work.
Meanwhile, there still remains the error you've brought to my
attention. The problem appears to be with the order in which the
poles/roots are sorted. The script residue.m calls mpoles.m, which
uses the internal sort() routine, to sort (reorder) the roots.
It appears that your octave reorders the roots differently than mine.
I've verified that the mpoles.m checked into the CVS is the same as
the one I'm using (however, there was a prior version I had posted
that would produce this problem).
I haven't yet had the time to mirror the entire CVS to my system. So
I am unable to isolate the problem, and be certain of my guess.
However, I think we can verify the hypothesis quickly. Using my FInk
install of Octave, I get the following
b = [1, 0, 1];
a = [1, 0, 18, 0, 81];
[r, p, k, e] = residue(b, a)
r =
0.00000 - 0.09259i
0.22222 + 0.00000i
0.00000 + 0.09259i
0.22222 - 0.00000i
p =
0.0000 + 3.0000i
0.0000 + 3.0000i
-0.0000 - 3.0000i
-0.0000 - 3.0000i
k = [](0x0)
e =
1
2
1
2
On your system, I suspect the poles are ordered as
p =
0.0000 - 3.0000i
0.0000 - 3.0000i
-0.0000 + 3.0000i
-0.0000 + 3.0000i
i.e. the -3i pair comes ahead of the +3i pair.
Can you verify?
If the problem is that simple, then I'll venture a guess that the
sort function has been modified since the 2.9.14 release ... or ...
your CPU is giving slightly different answers than mine. In either
event, it makes since to include both results in the test-script
(since they are each actually correct).
Assuming my assumptions are correct, I've modified the script. The
patch is below. Please let me know if the syntax for creating the
patch is not what you prefer.
diff -Naur residue_cvs.m residue.m
--- residue_cvs.m 2007-10-08 21:10:06.000000000 -0400
+++ residue.m 2007-10-08 22:06:07.000000000 -0400
@@ -122,31 +122,6 @@
## Created: June 1994
## Adapted-By: jwe
-%!test
-%! b = [1, 1, 1];
-%! a = [1, -5, 8, -4];
-%! [r, p, k, e] = residue (b, a);
-%! assert ((abs (r - [-2; 7; 3]) < 1e-5
-%! && abs (p - [2; 2; 1]) < 1e-7
-%! && isempty (k)
-%! && e == [1; 2; 1]));
-%! k = [1 0];
-%! [b, a] = residue (r, p, k);
-%! assert ((abs (b - [1, -5, 9, -3, 1]) < 1e-12
-%! && abs (a - [1, -5, 8, -4]) < 1e-12));
-
-%!test
-%! b = [1, 0, 1];
-%! a = [1, 0, 18, 0, 81];
-%! [r, p, k, e] = residue(b, a);
-%! assert ((abs (54*r - [-5i; 12; 5i; 12]) < 1e-6
-%! && abs (p - [3i; 3i; -3i; -3i]) < 1e-7
-%! && isempty (k)
-%! && e == [1; 2; 1; 2]));
-%! [br, ar] = residue (r, p, k);
-%! assert ((abs (br - b) < 1e-12
-%! && abs (ar - a) < 1e-12));
-
function [r, p, k, e] = residue (b, a, varargin)
if (nargin < 2 || nargin > 3)
@@ -323,12 +298,27 @@
endfunction
-%% test/octave.test/poly/residue-1.m
%!test
%! b = [1, 1, 1];
%! a = [1, -5, 8, -4];
%! [r, p, k, e] = residue (b, a);
-%! assert((abs (r - [-2; 7; 3]) < 1e-6
+%! assert ((abs (r - [-2; 7; 3]) < 1e-5
%! && abs (p - [2; 2; 1]) < 1e-7
%! && isempty (k)
%! && e == [1; 2; 1]));
+%! k = [1 0];
+%! [b, a] = residue (r, p, k);
+%! assert ((abs (b - [1, -5, 9, -3, 1]) < 1e-12
+%! && abs (a - [1, -5, 8, -4]) < 1e-12));
+
+%!test
+%! b = [1, 0, 1];
+%! a = [1, 0, 18, 0, 81];
+%! [r, p, k, e] = residue(b, a);
+%! assert ((abs (54*r - [-5i; 12; 5i; 12]) < 1e-6
+%! && abs (p - [3i; 3i; -3i; -3i]) < 1e-7
+%! && isempty (k)
+%! && e == [1; 2; 1; 2]));
+%! [br, ar] = residue (r, p, k);
+%! assert ((abs (br - b) < 1e-12
+%! && abs (ar - a) < 1e-12));
Ben
- Re: modified residues() for matlab compatibility, (continued)
- Re: modified residues() for matlab compatibility, John W. Eaton, 2007/10/05
- Re: modified residues() for matlab compatibility, Ben Abbott, 2007/10/05
- Re: modified residues() for matlab compatibility, John W. Eaton, 2007/10/05
- Re: modified residues() for matlab compatibility, Ben Abbott, 2007/10/05
- Re: modified residues() for matlab compatibility, John W. Eaton, 2007/10/06
- Re: modified residues() for matlab compatibility, John W. Eaton, 2007/10/06
- Re: modified residues() for matlab compatibility, Ben Abbott, 2007/10/06
- Re: modified residues() for matlab compatibility, Ben Abbott, 2007/10/06
- Message not available
- Re: modified residues() for matlab compatibility, Ben Abbott, 2007/10/06
- Re: modified residues() for matlab compatibility, John W. Eaton, 2007/10/08
- Re: modified residues() for matlab compatibility,
Ben Abbott <=
- Re: modified residues() for matlab compatibility, John W. Eaton, 2007/10/09
- Re: modified residues() for matlab compatibility, Ben Abbott, 2007/10/09
- Re: modified residues() for matlab compatibility, John W. Eaton, 2007/10/09