[Octave-bug-tracker] [bug #52919] comparison of complex values

[Michael Leitner]

Follow-up Comment #13, bug #52919 (project octave):

To provisionally sum up the different issues that have crept up in this
thread:

- the performance of sort() on complex values (compiled C++ is significantly
slower than trivial m-code solution)
- should the order on complex values be changed (either for principal reasons,
for least surprise, or for compatibility with Matlab, and if so, also for
sort/min/max)?
- should perhaps a warning be emitted, as comparing complex values is arguably
often unintentional?

And now a few thoughts on the comments over the weekend:

Re comment #12:
Yes, it is indeed slow. I think the difference is even larger than a factor
2.5. Consider this:

N=1e7;
r=floor(randn(N,1)+1i*randn(N,1))+.5+.5*1i;

tic;[~,j]=sortrows([abs(r) angle(r)]);rs1=r(j);toc
tic;rs2=sort(r);toc
isequal(rs1,rs2)

tic;[~,j]=unique([abs(r) angle(r)],"rows");u1=r(j);toc
tic;u2=unique(r);toc
isequal(u1,u2)


On my system, the difference is nearly a factor of four, and sortrows sorts
two times, not just once according to abs(r). This is just the code that is
also in sort_complex_test.m added in comment #6, using the fact that this is
actually the algorithm used by sortrows. With the rounding, also some absolute
values are equal, so that the decision on arg is relevant. Note that adding
0.5 to the numbers is necessary, as by conventional rounding you can end up
with -1-0i and -1+0i, which seemingly are treated differently (which probably
is also a small buglet). 

In C++, one could definitely go below the time taken by the sortrows-solution
by first sorting according to abs, and then make another pass and sort all
stretches with same abs by arg. This approach should take N* log(N)+N* log(N*
p), where p is the probability for two values having the same absolute value,
and which is potentially very small. On the other hand, sortrows takes 2N*
log(N) (disregarding the fact that also sort itself probably takes shortcuts
when many arguments are equal). But I think, the most sensible solution would
be just to change sort to do it in the sortrows-way for complex values. 

I support comment #11, with the proviso that I am not quite convinced that it
is a good idea to sort complex values (i.e., to rely on an order relation) in
any case. But yes, if there is an order relation, it should be consistent.

Re comment #10:
As regards the Matlab behaviour: probably it was not that bad, because I would
hope that the used sort algorithm is stable, that is, it keeps the input order
for values that compare as equal. According to the documentation, Octave uses
stable sort, which is also a prerequisite for sortrows.

As regards whether a real number is a special case of a complex number or not:
I do not know what Matlab does, but for Octave I would argue that it is.
Consider this:


a=[-1 1i];
iscomplex(a)
iscomplex(a(1))
iscomplex(a(2))


So without having told Octave to cast a(1) to real, it considers it to be
real. And of course, it uses the real ordering relations on it, which gives
the inconsistency that was my initial motivation for filing this bug report,
because the following two are different:


(a<1)(1)
a(1)<1


So if you want to hold up this concept of "Octave has two categories of
numbers", then you would need to have b=a(1);iscomplex(b) in the present
example evaluate to true. Personally, I would not do this.

So, now for the question which ordering to use on the complex numbers: For the
present (abs,arg)-sorting you can even add a real number to a complex number
to reverse the order, compare the results of (a<1) and (a+1<1+1). When you do
lexicographic (real,imag)-sorting, then this does not happen. Thus, with
lexicographic sorting you could consider the real numbers as special case of
complex numbers not only with respect to arithmetics, but also ordering, and
these issues that real-valued entries taken out of a complex vector compare
differently than the element-wise operation on the whole vector would not
appear.

The pseudocode in the comment actually just implements lexicographic ordering,
because if neither areal!=breal and aimag!=bimag, then of course a==b. I think
there is no way to "have an ordering that would be consistent with the real
numbers as the imaginary part of the number was shrunk to zero" but still sort
primarily on the absolute value, because this would always make -2+1i* eps
greater than 1+1i* eps.

Re comment #9:
Actually, when you want to efficiently sort complex numbers, you should rather
do it not in C++ (at least not in the present C++), as I have shown. I have no
idea at all what rlocus does and what typical input arguments are, but looking
into its code I would be surprised if the performance of the internal sort()
is actually deciding here. Are you sure that it is?

    _______________________________________________________

Reply to this item at:

  <http://savannah.gnu.org/bugs/?52919>

_______________________________________________
  Message sent via/by Savannah
  http://savannah.gnu.org/