[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: modified residues() for matlab compatibility
From: 
John W. Eaton 
Subject: 
Re: modified residues() for matlab compatibility 
Date: 
Tue, 09 Oct 2007 14:49:20 0400 
On 8Oct2007, Ben Abbott wrote:
 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 testscripts, I'd placed my testing between
 the docstrip 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?
Here is what I see with the current sources:
octave:1> b = [1, 0, 1];
octave:2> a = [1, 0, 18, 0, 81];
octave:3> [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
 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
I don't see any changes that should affect this result and I get the
sam result if I use the new residue with Octave 2.9.14.
 ... or ...
 your CPU is giving slightly different answers than mine. In either
 event, it makes since to include both results in the testscript
 (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.
This patch did not apply for me. Please use
cvs diff u
in the scripts/polynomial directory where you have Octave checked out.
jwe
 Re: modified residues() for matlab compatibility, (continued)
 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, 2007/10/08
 Re: modified residues() for matlab compatibility,
John W. Eaton <=
 Re: modified residues() for matlab compatibility, Ben Abbott, 2007/10/09
 Re: modified residues() for matlab compatibility, John W. Eaton, 2007/10/09