Re: isreal and iscomplex

[Daniel J Sebald]

On 09/10/2012 03:28 PM, Rik wrote:
> 9/10/12
>
> Dan,
>
> I think part of the confusion is that you are encountering two sets of
> functionality:
>
> 1) isreal function
> 2) automatic narrowing of complex->real
>
> The isreal function, as you now know, does not check that the imaginary
> part of an object is zero.  Rather it checks that the storage class used
> for the object is real (8 bytes) versus complex (16 bytes).  It is a hard
> requirement to keep isreal compatible with the Matlab function of the same
> name so we can't re-define it.  The other functions, iscomplex and isimag,
> are Octave's alone and could be re-defined.  As Michael has said, there
> would need to be a very clear plan laid out so that it doesn't compromise
> older scripts, etc.
>
> The second functionality that you are encountering is the automatic
> narrowing that Octave performs on complex numbers.  Complex numbers take
> twice as much memory to store and there is a significant incentive to
> reduce results to real numbers if possible.  Octave does this at many
> places: when you type a number on the command line, when you define a
> matrix, when you extract a value or range using indexing, etc.  Thus,
>
> isreal (1+0i) =>  isreal (1) =>  1
> isreal ([1+0i, 1]) =>  isreal ([1, 1]) =>  1
> x = [1, 1i, 1+1i];
> isreal (x) =>  false because x is stored as a complex object
> isreal (x(1)) =>  isreal (1) =>  true
>
> The complex() function is guaranteed to return a complex result.  Therefore,
>
> isreal (complex(1,0)) =>  1 because the object is complex EVEN though
>                             the imaginary part is zero.
>
> However, it is simple enough to re-engage narrowing by using the matrix
> constructor operator '['.
>
> isreal ([complex(1,0)]) =>  isreal ([1 + 0i]) =>  isreal ([1]) =>  1

OK, I'm starting to get it now.  Thanks.


> Incidentally, I did a little benchmarking and automatic narrowing is ~25%
> faster than using an equivalent vectorized operator.
>
> Example code:
> x = complex (ones (1e4), zeros (1e4));
>
> tic; xreal = isreal ([x]); toc
> Elapsed time is 1.78048 seconds.
>
> tic; xreal = all (imag (x) == 0); toc
> Elapsed time is 2.42 seconds.

That's not what I'm getting...

octave:1> x = complex (ones (1e4), zeros (1e4));
octave:2>
octave:2> cpustart = cputime();
octave:3> xisreal = isreal(x)
xisreal = 0
octave:4> cputime() - cpustart
ans =  0.0029990
octave:5>
octave:5> cpustart = cputime();
octave:6> xisreal = isreal([x])
xisreal =  1
octave:7> cputime() - cpustart
ans =  1.5658
octave:8>
octave:8> cpustart = cputime();
octave:9> xisreal = isreal(x(:))
xisreal =  1
octave:10> cputime() - cpustart
ans =  1.5678
octave:11>
octave:11> cpustart = cputime();
octave:12> xisreal = all (all (imag (x) == 0))
xisreal =  1
octave:13> cputime() - cpustart
ans =  1.4528

What brand of computer do you have there Rik?


> In your patch for fftfilt you might want to take advantage of that fact.
> Also, if you are already using indexing then the narrowing will happen
> automatically.  I could test whether column 2 of x is truly real (imaginary
> part == 0) with
>
> isreal (x(:,2))
>
> which reads more naturally then
> all (imag (x(:,2)) == 0)

I can straighten that out now given I understand how isreal() works.


> You might want to check Matlab's own documentation for isreal
> (http://www.mathworks.com/help/techdoc/ref/isreal.html).  They do a good
> job in the introduction and the Tips section of explaining all these corner
> cases.
>
> Octave could certainly use more in the documentation on this behavior since
> it can be confusing.

Yes, definitely more documentation than just the two sentences currently 
given.  Having sorted through all this I can now summarize as follows:

...

The function isreal(A) tests that the STORAGE CLASS of variable A uses 
only a real component.  However, storage class IS NOT PRESERVED through 
any mathematical operations.  That is, reductions in memory space are 
attempted with all matrix computations.  Therefore, do not use isreal(A) 
to test if all imaginary components of a complex matrix are zero. 
Instead, use isreal(A(:)) or all(all(imag(A)==0)).

...

Well, now let's return to the notion of having an "iscomplex()" and 
"isimag()".  Is there some utility to that?  I conclude "Maybe".  Given 
that "isreal" is cast in this storage class paradigm, the only time it 
would make sense to have an "isimag()" is if there were a similar 
reduction to 8-byte matrices if all real parts of the complex matrix are 
zero.  (To have isimag() mean otherwise would REALLY confuse matters.) 
Currently, I'm guessing Octave isn't set up that way, i.e., to save 
space after a mathematical operation by using just an imaginary 
component if all real components are zero.

One little non-symmetry for a real-only/imaginary-only paradigm would be 
the case where the complex matrix operation result is all 0 + 0i, in 
which case the storage class would be "real" vector and all zeros.  With 
that behavior, "isreal()" would still behave in exactly the same way as 
it currently does.  However, note that "iscomplex()" would not be the 
complement of "isreal()" under such a scenario because if the math 
result had all real components zero and were stored in the 8-byte 
imaginary component only then iscomplex() would be false.  (And note I 
avoided the word "complex" in my documentation summary above.)

Personally, I'd prefer "iscomplex()" be deprecated in the short term 
since it currently is only the complement of "isreal()".  That would 
leave open its use in a real/imaginary/complex storage class 
configuration, which is of very low priority right now.  It would also 
mean we have but one hunk of tricky function documentation to worry 
about on the matter.

Dan