bugGNU Octave - Bugs: bug #50124, cplxpair produces incorrect result...

 
 

bug #50124: cplxpair produces incorrect result for second column of complex pairs

Submitter:  Bill Lash <welash>
Submitted:  Mon 23 Jan 2017 12:07:27 AM UTC
   
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Incorrect Result
Status:  Fixed Assigned to:  None
Originator Name:  William Lash Open/Closed:  * Closed
Release:  * 4.0.3 Operating System:  * GNU/Linux
Fixed Release:  None Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

Jump to the original submission

Wed 30 Aug 2017 12:37:24 AM UTC, comment #25: 

I tested the script with the changes that Rik made, and all of my test cases pass.  I think it can be closed.  Thanks!

Bill Lash <welash>
Tue 29 Aug 2017 10:46:41 PM UTC, comment #24: 

I made several changes to update the documentation and follow Octave coding conventions here (http://hg.savannah.gnu.org/hgweb/octave/rev/a315ac23dc6c).

The only thing of note was this change


-  midx=idx+[0:si(2)-1]*si(1);
+  midx = idx + rows (idx) * ones (rows (idx), 1) * [0:columns(idx)-1];


This was done to be analogous to


  z = z(idx + n * ones (n, 1) * [0:m-1]);


If we don't need the call to ones(), because we are relying on implicit broadcasting, then it can be removed from this line as well.

Rik <rik5>
Group administrator
Tue 29 Aug 2017 09:08:37 PM UTC, comment #23: 

I pushed the change here:

http://hg.savannah.gnu.org/hgweb/octave/rev/b4cf41d173dd

Can this report be closed?

John W. Eaton <jwe>
Group administrator
Mon 17 Jul 2017 05:50:55 PM UTC, comment #22: 

Nothing to be done.  It just needs someone to do the final review and push.  But it's the middle of summer and even Octave volunteers take time off.  Also, there has been some really deep re-structuring of the Octave code base around using C++11 features and that has sucked up a lot of developer time.

Rik <rik5>
Group administrator
Sun 16 Jul 2017 01:31:12 AM UTC, comment #21: 

I happened to think about this the other day, and notice that the patch has not yet been applied.  Is there anything I need to do?  I think the latest version obeys coding standards and has tests.

Bill Lash <welash>
Sun 29 Jan 2017 05:38:29 AM UTC, comment #20: 

Thanks for matlab-izing the test script,  I always forget that octave has a superset of the matlab syntax.

So as far as I can tell, these are the issues where matlab and octave disagree:

  1.  The issue reported originally, when passed an array with columns of complex pairs, the results for all but the first column are incorrect


  1.  When passsing a tolerance value of 1 or greater, matlab returns an error


  1.  Values that cplxpair determines to be real (i.e. have an imaginary part less than tol * abs(z)) are forced to be real


  1.  Complex pairs are forced to be true complex conjugates.  It appears that Matlab chooses the member of the pair that has a positive imaginary part, and conjugates it to get the member of the pair with the negative part.  The pairs are sorted with the member that has negative real part listed first.


  1.  The way that matlab seems to be comparing values to determine if they are complex pairs is by comparing the conjugate of each value to the other values, and finding the one with the minimum magnitude of the difference.  If this minimum difference is less than tol * abs(z) then it is the pair. 


I tried my hand at adressing these issues, and I think I have version that works.  I made the changes against the latest version in mercurial.  I am attaching a patch (new_cplxpair.patch) as well as the results of running your version of the script (latest.out).  There are some formatting differences (matlab gives one more digit than octave in format long, and in some cases octave prints 1 or 0 instead of 1.000000000000000), but it seems to match up other than that.

I did have to change some of the tests, due to forcing of complex pairs to be conjugates causing test failures.

(file #39599, file #39600)

Bill Lash <welash>
Fri 27 Jan 2017 02:25:39 PM UTC, comment #19: 

for future reference, had to cleanup the script to run in Matlab:

no printf in matlab.

printf --> fprintf

" --> '

endfor --> end

end_try_catch --> end

"Undefined function 'printf' for input arguments of type 'char'"

can't do abs(z)(1)

   absz = abs(z)
   absz(1)


anyway, after that it ran and I don't think I changed any intended function:


>> tcplxpair
Tests to determine how the tolerance parameter is handled for determining real numbers
How is the default tolerance used?
Was real until 1 + 100*eps*i
abs(z) = 1
Do we get the same results by passing 100*eps for the tolerance?
Was real until 1 + 100*eps*i
abs(z) = 1
Does doubling the tolerance work as expected?
Was real until 1 + 200*eps*i
abs(z) = 1
Tests with default tolerance
Determine how real numbers are found
Was real until 1 + 100*eps*i
abs(z) = 1
If the real part is doubled how does it compare
Was real until 2 + 200*eps*i
abs(z) = 2
What if there is another real value of large magnitude
Was real until 2 + 200*eps*i
abs(z) = 2
Tests for complex pairs
What if there is error in the real part, but imaginary parts are equal
Pairing was successful until (1 + 141*eps) + 1i, 1 - 1i
abs(z) = 1.41421
if the imaginary part is doubled what happens
Pairing was successful until (1 + 223*eps) + 2i, 1 - 2i
abs(z) = 2.23607
how are errors in the imaginary part handled
Pairing was successful until 1 + (141*eps*1i + 1i), 1 - 1i
abs(z) = 1.41421
scale the real part
Pairing was successful until 2 + (223*eps*1i + 1i), 2 - 1i
abs(z) = 2.23607
scale the imaginary part
Pairing was successful until 1 + (222*eps*1i + 2i), 1 - 2i
abs(z) = 2.23607
Add a larger pair of cojugates to make sure it is being done per pair
Pairing was successful until 1 + (222*eps*1i + 2i), 1 - 2i, 400 + 100i, 400 - 100i
abs(z) = 2.23607
sort order with a real, and error in the imaginary part
z =
  1.000000000000000 + 0.000000000000000i
  1.000000000000000 + 1.000000000000018i
  1.000000000000000 - 1.000000000000000i
ans =
  1.000000000000000 - 1.000000000000018i
  1.000000000000000 + 1.000000000000018i
  1.000000000000000 + 0.000000000000000i
z =
  1.000000000000000 + 0.000000000000000i
  1.000000000000000 + 1.000000000000000i
  1.000000000000000 - 1.000000000000018i
ans =
  1.000000000000000 - 1.000000000000000i
  1.000000000000000 + 1.000000000000000i
  1.000000000000000 + 0.000000000000000i
does the input order matter for the real part
z =
  1.000000000000000 + 1.000000000000018i
  1.000000000000000 - 1.000000000000000i
  1.000000000000000 + 0.000000000000000i
ans =
  1.000000000000000 - 1.000000000000018i
  1.000000000000000 + 1.000000000000018i
  1.000000000000000 + 0.000000000000000i
z =
  1.000000000000000 + 1.000000000000000i
  1.000000000000000 - 0.999999999999982i
  1.000000000000000 + 0.000000000000000i
ans =
  1.000000000000000 - 1.000000000000000i
  1.000000000000000 + 1.000000000000000i
  1.000000000000000 + 0.000000000000000i
does the input order matter for the cplx pairs
z =
  1.000000000000018 + 1.000000000000000i
  1.000000000000000 - 1.000000000000000i
ans =
  1.000000000000018 - 1.000000000000000i
  1.000000000000018 + 1.000000000000000i
z =
  1.000000000000000 - 1.000000000000000i
  1.000000000000018 + 1.000000000000000i
ans =
  1.000000000000018 - 1.000000000000000i
  1.000000000000018 + 1.000000000000000i
If value is real, with imaginary component less than tolerance, is result real?z =
  1.000000000000000 + 0.000000000000018i
ans =
     1
ans =
  logical
   1


i attached the matlab compatible version as tcplxpair_matlab.m

(file #39570)

Nicholas Jankowski <nrjank>
Group Member
Fri 27 Jan 2017 06:33:42 AM UTC, comment #18: 

Should there also be a dependency on bug #47865 https://savannah.gnu.org/bugs/?47865 that was mentioned earlier, which seems to be trying to address the tolerance issues.

If matlab supports try/catch I came up with a test script the ranges over a set of "error" values for real numbers and complex pairs to try to determine how the tolerance argument is used in each case. 

It also does a few tests to determine how the results are sorted and which value of the complex pair is used when forcing the complex pairs to be equal.  I am attaching the script, called tcplxpair, and the results of running it as:

octave --no-gui tcplxpair > output

I ran it with 4.0.3, and I also downloaded the current version of the cplxpair.m file from mercurial (which seems to also be in octave 4.2) and got the results from that.

If someone could try it with matlab and share the results, I'll try to analyze the differences.

(file #39566, file #39567, file #39568)

Bill Lash <welash>
Fri 27 Jan 2017 05:18:48 AM UTC, comment #17: 

I added a dependency on bug #45810.  That was the bug that motivated the last serious change to the tolerance calculation.

Rik <rik5>
Group administrator
Fri 27 Jan 2017 05:01:53 AM UTC, comment #16: 

So it now seems that there are multiple incompatibilities.  First, the fact that the second column of values is not sorted correctly.  Second, for Matlab compatibility, results should be exact conjugates.  This is documented by The Mathworks:


The complex conjugate pairs are forced to be exact complex conjugates


Third, the tolerance needs to be in the range 0 <= tol <= 1 and probably the algorithm in cplxpair.m updated as well.


Rik <rik5>
Group administrator
Fri 27 Jan 2017 04:30:40 AM UTC, comment #15: 


>> z = [1 + 0.5i; 1; 1-0.49i]

z =

  1.000000000000000 + 0.500000000000000i
  1.000000000000000 + 0.000000000000000i
  1.000000000000000 - 0.490000000000000i

>> cplxpair (z, 0.1)

ans =

  1.000000000000000 - 0.500000000000000i
  1.000000000000000 + 0.500000000000000i
  1.000000000000000 + 0.000000000000000i




>> zz=exp (2i*pi*[0:4]'/5)

zz =

  1.000000000000000 + 0.000000000000000i
  0.309016994374947 + 0.951056516295154i
 -0.809016994374947 + 0.587785252292473i
 -0.809016994374947 - 0.587785252292473i
  0.309016994374947 - 0.951056516295154i

>> cplxpair(zz) - cplxpair(conj(zz))

ans =

   1.0e-15 *

  0.111022302462516 - 0.222044604925031i
  0.111022302462516 + 0.222044604925031i
  0.222044604925031 + 0.111022302462516i
  0.222044604925031 - 0.111022302462516i
  0.000000000000000 + 0.000000000000000i


Nicholas Jankowski <nrjank>
Group Member
Fri 27 Jan 2017 04:22:04 AM UTC, comment #14: 

Can you also try:


z = [1 + 0.5i; 1; 1-0.49i]
cplxpair (z, 0.1)


From your results it looks like tolerance needs to be between 0 and 1, and is not multiplied by eps internally (since (1-eps)*eps*abs(z) would not be greater 0.01)

Also, for the FIXME, can you try:


zz=exp (2i*pi*[0:4]'/5);
cplxpair(zz) - cplxpair(conj(zz))



Also, some of this is getting away from the original issue that I reported.  Should I create another bug report?

Bill Lash <welash>
Fri 27 Jan 2017 03:45:33 AM UTC, comment #13: 

Re: the FIXME question in Comment #7

matlab 2016b:


>> exp (2i*pi*[0:4]'/5)

ans =

  1.000000000000000 + 0.000000000000000i
  0.309016994374947 + 0.951056516295154i
 -0.809016994374947 + 0.587785252292473i
 -0.809016994374947 - 0.587785252292473i
  0.309016994374947 - 0.951056516295154i

>> cplxpair(ans)

ans =

 -0.809016994374947 - 0.587785252292473i
 -0.809016994374947 + 0.587785252292473i
  0.309016994374947 - 0.951056516295154i
  0.309016994374947 + 0.951056516295154i
  1.000000000000000 + 0.000000000000000i


Nicholas Jankowski <nrjank>
Group Member
Fri 27 Jan 2017 03:42:14 AM UTC, comment #12: 

Re: Comment #7

Matlab 2016b:

>> z = [1 + 0.5i; 1; 1-0.49i]
cplxpair (z, 1e14)

z =

   1.0000 + 0.5000i
   1.0000 + 0.0000i
   1.0000 - 0.4900i

Error using cplxpair (line 41)
Relative tolerance TOL must be a scalar such that 0<=TOL<1.

>> cplxpair (z, 1)
Error using cplxpair (line 41)
Relative tolerance TOL must be a scalar such that 0<=TOL<1.

>> cplxpair (z, 1-eps)

ans =

     1
     1
     1


Nicholas Jankowski <nrjank>
Group Member
Fri 27 Jan 2017 01:53:42 AM UTC, comment #11: 

is this a recent change in behavior for octave?  With 4.0.3 I see:


>> z = [1 + 0.5i; 1; 1-0.49i]
z =

   1.0000 + 0.5000i
   1.0000 + 0.0000i
   1.0000 - 0.4900i

>> cplxpair (z, 1e14)
ans =

   1.0000 + 0.5000i
   1.0000 + 0.0000i
   1.0000 - 0.4900i

>> cplxpair (z, 0.1)
ans =

   1.0000 - 0.4900i
   1.0000 + 0.5000i
   1.0000 + 0.0000i


So for 4.0.3 it looks like cplxpair is comparing to tol*abs(z).

Bill Lash <welash>
Thu 26 Jan 2017 11:10:45 PM UTC, comment #10: 

To address Mike's comment: if we don't know Matlab's exact tolerance calculation, then start with a low tolerance.  With a low tolerance, the two complex numbers won't pair and Matlab will emit an error.  Increase the tolerance until Matlab is able to pair the two complex numbers and then see whether the result is an exactly conjugated pair or the original two complex numbers paired.

To address Bill's comment: The calculated tolerance used by Octave is

tol * eps (abs (Z))


For Z == 1+0.5i,


abs (1+0.5i)
ans =  1.1180
eps (ans)
ans =    2.2204e-16
ans * 1e14
ans =  0.022204


The test metric is


abs (Z2 - conj (Z1))


or


abs (1-0.49i - conj (1+0.5i))
===
abs (1-0.49i - (1-0.5i))
===
abs (.01i)
===
.01


Since .01 < 0.022204, they will pair.

If you introduce an absolutely crazy tolerance, like 1e25, then the numbers will be treated as all real.  However, this is easily detected because the true real number (1+0i) will not be at the end of the list separated from the complex numbers at the start.  Instead, the list ordering will be the same as the input in Octave.

Rik <rik5>
Group administrator
Thu 26 Jan 2017 10:25:28 PM UTC, comment #9: 

Won't setting the tolerance that high mean that all of the data you have given will be treated as real?  I would think using a tolerance of 0.1 would be more appropriate.  So


z = [1 + 0.5i; 1; 1-0.49i]
cplxpair (z, 0.1)



Bill Lash <welash>
Thu 26 Jan 2017 10:21:52 PM UTC, comment #8: 

We also have the still open question about how that tolerance argument is actually interpreted, since Matlab's docs are very confusing about it (bug #47865).

Mike Miller <mtmiller>
Group Member
Thu 26 Jan 2017 10:14:41 PM UTC, comment #7: 

The function cplxpair is defined in an m-file.  If you look at the start of the code there is a series of FIXME notes.


## FIXME: subsort returned pairs by imaginary magnitude
## FIXME: Why doesn't exp (2i*pi*[0:4]'/5) produce exact conjugates?  Does
## FIXME: it in Matlab?  The reason is that complex pairs are supposed
## FIXME: to be exact conjugates, and not rely on a tolerance test.


This seems relevant.  Also, with a high enough tolerance you can force Octave to match pairs which are fairly far apart.  Octave does pair them, but returns each half of the pair exactly as it appeared in the input.  If Matlab is choosing to keep only one half of a pair and conjugate it to get the second half, then this test will immediately show it.


z = [1 + 0.5i; 1; 1-0.49i]
cplxpair (z, 1e14)


The results for Octave are


octave:21> z = [1 + 0.5i; 1; 1-0.49i]
z =

   1.0000 + 0.5000i
   1.0000 + 0.0000i
   1.0000 - 0.4900i

octave:22> cplxpair (z, 1e14)
ans =

   1.0000 - 0.4900i
   1.0000 + 0.5000i
   1.0000 + 0.0000i


@Nicholas: Can you run the test code in Matlab to see what happens?


Rik <rik5>
Group administrator
Mon 23 Jan 2017 04:07:58 AM UTC, comment #6: 

Thanks, so it looks like Matlab is choosing one of the two values of the set of complex pairs, and conjugating it to get the other value.

I wonder if that is a big enough difference to file a matlab compatability bug against.

Bill Lash <welash>
Mon 23 Jan 2017 03:55:58 AM UTC, comment #5: 

Matlab 2016b output:


>> format long
>> z = exp (2i*pi*[4; 3; 5; 2; 6; 1; 0]/7);
>> y = [ -1-1i; -1+1i;-3; -2; 1; 2; 3];
>> [z,y]
ans =
 -0.900968867902419 - 0.433883739117558i -1.000000000000000 - 1.000000000000000i
 -0.900968867902419 + 0.433883739117558i -1.000000000000000 + 1.000000000000000i
 -0.222520933956315 - 0.974927912181824i -3.000000000000000 + 0.000000000000000i
 -0.222520933956314 + 0.974927912181824i -2.000000000000000 + 0.000000000000000i
  0.623489801858733 - 0.781831482468030i  1.000000000000000 + 0.000000000000000i
  0.623489801858734 + 0.781831482468030i  2.000000000000000 + 0.000000000000000i
  1.000000000000000 + 0.000000000000000i  3.000000000000000 + 0.000000000000000i
>> cplxpair([z,y])
ans =
 -0.900968867902419 - 0.433883739117558i -1.000000000000000 - 1.000000000000000i
 -0.900968867902419 + 0.433883739117558i -1.000000000000000 + 1.000000000000000i
 -0.222520933956314 - 0.974927912181824i -3.000000000000000 + 0.000000000000000i
 -0.222520933956314 + 0.974927912181824i -2.000000000000000 + 0.000000000000000i
  0.623489801858734 - 0.781831482468030i  1.000000000000000 + 0.000000000000000i
  0.623489801858734 + 0.781831482468030i  2.000000000000000 + 0.000000000000000i
  1.000000000000000 + 0.000000000000000i  3.000000000000000 + 0.000000000000000i
>> cplxpair([z,y]) - [z,y]
ans =
   1.0e-15 *
  0.111022302462516 - 0.222044604925031i  0.000000000000000 + 0.000000000000000i
  0.000000000000000 + 0.000000000000000i  0.000000000000000 + 0.000000000000000i
  0.249800180540660 + 0.000000000000000i  0.000000000000000 + 0.000000000000000i
  0.000000000000000 + 0.000000000000000i  0.000000000000000 + 0.000000000000000i
  0.222044604925031 + 0.111022302462516i  0.000000000000000 + 0.000000000000000i
  0.000000000000000 + 0.000000000000000i  0.000000000000000 + 0.000000000000000i
  0.000000000000000 + 0.000000000000000i  0.000000000000000 + 0.000000000000000i


Nicholas Jankowski <nrjank>
Group Member
Mon 23 Jan 2017 03:52:06 AM UTC, comment #4: 

Thanks, that is more what I expected.  It is strange that matlab shows a slight difference between the before and after for the first set of complex pairs.  Would it be possible to do that again with "format long"?

It seems like it is possible that matlab is forcing the complex pairs to be conjugates of each other.  With my proposed patch, and adding a format long at the front of the script, I get:


>> format long
>> z = exp (2i*pi*[4; 3; 5; 2; 6; 1; 0]/7);
>> y = [ -1-1i; -1+1i;-3; -2; 1; 2; 3];
>> [z,y]
ans =

  -0.900968867902419 - 0.433883739117558i  -1.000000000000000 - 1.000000000000000i
  -0.900968867902419 + 0.433883739117558i  -1.000000000000000 + 1.000000000000000i
  -0.222520933956315 - 0.974927912181824i  -3.000000000000000 + 0.000000000000000i
  -0.222520933956314 + 0.974927912181824i  -2.000000000000000 + 0.000000000000000i
   0.623489801858733 - 0.781831482468030i   1.000000000000000 + 0.000000000000000i
   0.623489801858734 + 0.781831482468030i   2.000000000000000 + 0.000000000000000i
   1.000000000000000 + 0.000000000000000i   3.000000000000000 + 0.000000000000000i

>> cplxpair([z,y])
ans =

  -0.900968867902419 - 0.433883739117558i  -1.000000000000000 - 1.000000000000000i
  -0.900968867902419 + 0.433883739117558i  -1.000000000000000 + 1.000000000000000i
  -0.222520933956315 - 0.974927912181824i  -3.000000000000000 + 0.000000000000000i
  -0.222520933956314 + 0.974927912181824i  -2.000000000000000 + 0.000000000000000i
   0.623489801858733 - 0.781831482468030i   1.000000000000000 + 0.000000000000000i
   0.623489801858734 + 0.781831482468030i   2.000000000000000 + 0.000000000000000i
   1.000000000000000 + 0.000000000000000i   3.000000000000000 + 0.000000000000000i

>> cplxpair([z,y]) - [z,y]
ans =

   0   0
   0   0
   0   0
   0   0
   0   0
   0   0
   0   0


Bill Lash <welash>
Mon 23 Jan 2017 03:41:47 AM UTC, comment #3: 

matlab 2016b output for the test script:


>> z = exp (2i*pi*[4; 3; 5; 2; 6; 1; 0]/7);
>> y = [ -1-1i; -1+1i;-3; -2; 1; 2; 3];
>> [z,y]
ans =
  -0.9010 - 0.4339i  -1.0000 - 1.0000i
  -0.9010 + 0.4339i  -1.0000 + 1.0000i
  -0.2225 - 0.9749i  -3.0000 + 0.0000i
  -0.2225 + 0.9749i  -2.0000 + 0.0000i
   0.6235 - 0.7818i   1.0000 + 0.0000i
   0.6235 + 0.7818i   2.0000 + 0.0000i
   1.0000 + 0.0000i   3.0000 + 0.0000i
>> cplxpair([z,y])
ans =
  -0.9010 - 0.4339i  -1.0000 - 1.0000i
  -0.9010 + 0.4339i  -1.0000 + 1.0000i
  -0.2225 - 0.9749i  -3.0000 + 0.0000i
  -0.2225 + 0.9749i  -2.0000 + 0.0000i
   0.6235 - 0.7818i   1.0000 + 0.0000i
   0.6235 + 0.7818i   2.0000 + 0.0000i
   1.0000 + 0.0000i   3.0000 + 0.0000i
>> cplxpair([z,y]) - [z,y]
ans =
   1.0e-15 *
   0.1110 - 0.2220i   0.0000 + 0.0000i
   0.0000 + 0.0000i   0.0000 + 0.0000i
   0.2498 + 0.0000i   0.0000 + 0.0000i
   0.0000 + 0.0000i   0.0000 + 0.0000i
   0.2220 + 0.1110i   0.0000 + 0.0000i
   0.0000 + 0.0000i   0.0000 + 0.0000i
   0.0000 + 0.0000i   0.0000 + 0.0000i


Nicholas Jankowski <nrjank>
Group Member
Mon 23 Jan 2017 12:43:07 AM UTC, comment #2: 

I think the patch below fixes the issue.  I am curious if there is a better way to rearrange the z matrix than the way I did it.  Basically idx contains columns of numbers between 1 and 7 in my test case, and the intent is to rearrange each column of z according to the corresponding column of idx.  I created midx which creates an array of pairs of coordinates.

Comments are welcome.

(file #39534)

Bill Lash <welash>
Mon 23 Jan 2017 12:32:11 AM UTC, comment #1: 

oops, the cut and paste lost something in the markup.  The test script is actually:


z = exp (2i*pi*[4; 3; 5; 2; 6; 1; 0]/7);
y = [ -1-1i; -1+1i;-3; -2; 1; 2; 3];
[z,y]
cplxpair([z,y])
cplxpair([z,y]) - [z,y]


Bill Lash <welash>
Mon 23 Jan 2017 12:07:27 AM UTC, original submission:  

If I understand the documentation for cplxpair correctly, if you pass in a matrix with 2 columns, the function should sort the pairs in each column (unless you pass the DIM argument to tell it to sort along the specified dimension).

I am seeing that the second column is filled in with values from the first column.  Here is my test script:

z = exp (2i*pi*[4; 3; 5; 2; 6; 1; 0]/7);
y = [ -1-1i; -1+1i;-3; -2; 1; 2; 3];
[z,y]
cplxpair([z,y])
cplxpair([z,y]) - [z,y]


I get:

>> z = exp (2i*pi*[4; 3; 5; 2; 6; 1; 0]/7);
>> y = [ -1-1i; -1+1i;-3; -2; 1; 2; 3];
>> [z,y]

ans =

  -0.90097 - 0.43388i  -1.00000 - 1.00000i
  -0.90097 + 0.43388i  -1.00000 + 1.00000i
  -0.22252 - 0.97493i  -3.00000 + 0.00000i
  -0.22252 + 0.97493i  -2.00000 + 0.00000i
   0.62349 - 0.78183i   1.00000 + 0.00000i
   0.62349 + 0.78183i   2.00000 + 0.00000i
   1.00000 + 0.00000i   3.00000 + 0.00000i

>> cplxpair([z,y])

ans =

  -0.90097 - 0.43388i  -0.22252 - 0.97493i
  -0.90097 + 0.43388i  -0.22252 + 0.97493i
  -0.22252 - 0.97493i  -0.90097 - 0.43388i
  -0.22252 + 0.97493i  -0.90097 + 0.43388i
   0.62349 - 0.78183i   0.62349 - 0.78183i
   0.62349 + 0.78183i   0.62349 + 0.78183i
   1.00000 + 0.00000i   1.00000 + 0.00000i

>> cplxpair([z,y]) - [z,y]

ans =

   0.00000 + 0.00000i   0.77748 + 0.02507i
   0.00000 + 0.00000i   0.77748 - 0.02507i
   0.00000 + 0.00000i   2.09903 - 0.43388i
   0.00000 + 0.00000i   1.09903 + 0.43388i
   0.00000 + 0.00000i  -0.37651 - 0.78183i
   0.00000 + 0.00000i  -1.37651 + 0.78183i
   0.00000 + 0.00000i  -2.00000 + 0.00000i

>>


looking at the code in cplxpair.m, I am suspicious of a couple of things.  In:

  [idxi, idxj] = find (abs (imag (z)) ./ (abs (z) + realmin (cls)) < tol);
  q = sparse (idxi, idxj, 1, n, m);
  nr = sum (q, 1);
  [q, idx] = sort (q, 1);
  z = z(idx);
  y = z;

if I take out the semicolons, I see

idx =

   1   3
   2   4
   3   1
   4   2
   5   5
   6   6
   7   7

and then

z(idx) is

z =

  -0.90097 - 0.43388i  -0.22252 - 0.97493i
  -0.90097 + 0.43388i  -0.22252 + 0.97493i
  -0.22252 - 0.97493i  -0.90097 - 0.43388i
  -0.22252 + 0.97493i  -0.90097 + 0.43388i
   0.62349 - 0.78183i   0.62349 - 0.78183i
   0.62349 + 0.78183i   0.62349 + 0.78183i
   1.00000 + 0.00000i   1.00000 + 0.00000i

I think the indexing for the second column needs  to be relative to the second column.

Additionally in the loop that does the pairing, the variable j is not used in indexing the y or z array.


Bill Lash <welash>

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:
   
   
Comment:
   

Attached Files
file #39599:  new_cplxpair.patch added by welash (3KiB - text/x-patch - new version of patch to address all issues, and test results using this patch)
file #39600:  latest.out added by welash (3KiB - application/octet-stream - new version of patch to address all issues, and test results using this patch)
file #39570:  tcplxpair_matlab.m added by nrjank (4KiB - application/octet-stream - matlab compatible version of testscript)
file #39566:  tcplxpair added by welash (4KiB - application/octet-stream)
file #39567:  results4.2 added by welash (3KiB - application/octet-stream)
file #39568:  results_4.0.3 added by welash (3KiB - application/octet-stream)
file #39534:  cplxpair.patch added by welash (2KiB - text/x-patch)

 

Digest:
   bug dependencies.

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by jwe (Posted a comment)
  • -email is unavailable- added by rik5 (Posted a comment)
  • -email is unavailable- added by nrjank (Posted a comment)
  • -email is unavailable- added by welash (Submitted the item)
  •  

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

    Only group members can vote.

     

    Follow 14 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2019-03-03 mtmiller Carbon-CopyRemoved 80942 -
    2017-08-30 rik5 StatusReady For Test Fixed
        Open/ClosedOpen Closed
    2017-08-29 jwe StatusPatch Submitted Ready For Test
    2017-01-29 welash Attached File- Added new_cplxpair.patch, #39599
        Attached File- Added latest.out, #39600
    2017-01-27 nrjank Attached File- Added tcplxpair_matlab.m, #39570
    2017-01-27 siko1056 Carbon-CopyRemoved 91577 -
    2017-01-27 welash Attached File- Added tcplxpair, #39566
        Attached File- Added results4.2, #39567
        Attached File- Added results_4.0.3, #39568
    2017-01-27 rik5 Dependencies- Depends on bugs #45810
    2017-01-23 siko1056 StatusNone Patch Submitted
    2017-01-23 welash Attached File- Added cplxpair.patch, #39534

    Back to the top

    Powered by Savane 3.13-f8d8.
    Corresponding source code