January 02, 2014
On Wednesday, 1 January 2014 at 12:29:35 UTC, Stewart Gordon
wrote:
> On 19/11/2013 02:03, Andrei Alexandrescu wrote:> On 11/18/13 5:44 PM, Craig Dillabaugh wrote:
> <snip>
>>> Is there any reason why complex numbers in D's standard lib must
>>> be of non-integral types?  I believe in C++ the type is optimized
>>> for floating point values, but allows other types.
>>
>> The simple reason is we couldn't find appropriate applications.
> <snip>
>
> I don't understand.  At the moment Complex appears to me to be type-agnostic - as long as a type supports the standard arithmetic operators and assignment of the value 0 to it, it will work.  The only thing preventing it from working at the moment is this line
>
>     struct Complex(T)  if (isFloatingPoint!T)
>
> So why do you need an appropriate application in order not to have this arbitrary restriction?  Or have I missed something?
>
> It isn't just integer types.  Somebody might want to use complex with a library (fixed-point, arbitrary precision, decimal, etc.) numeric type.
>  Fractal generators, for example, are likely to use this a lot.
>
> Or even more exotically, use Complex!(Complex!real) to implement hypercomplex numbers.
>
> Stewart.

Since I originally proposed this I should chime in again. My
use-case for non floating point, complex numbers is radar image
processing, where the radar intensity/phase are stored as complex
numbers which are quantized as 16-bit integer values.

Andrei suggested I come up with a proposal for non-FP complex,
but reading this thread, and my experience working with integer
valued complex values in C++ has now reversed my opinion.  I
think the current D requirement is good.

Just a small example, the norm() method calculates the squared
magnitude of the complex number (eg. you get a single real
value).  This is an operation that we use extensively in our
work.  This created problems with generating lots of NaN values
due to integer overflows, that were somewhat tricky to find.

I can imagine numerous issues like this popping up if using
integral values with complex numbers.
January 02, 2014
On 2014-01-01 12:29, Stewart Gordon wrote:
> Or even more exotically, use Complex!(Complex!real) to implement
> hypercomplex numbers.

For a more generic solution to this, see Cayley-Dickson construction[1] and my implementation of such:

https://github.com/Biotronic/Collectanea/blob/master/biotronic/CayleyDickson2.d

It's nowhere near finished, has some bugs, and I haven't worked on it for at least a year, but it's an interesting way to create complex, hypercomplex, dual, and split-complex numbers (and combinations thereof).

[1]: http://en.wikipedia.org/wiki/Cayley-Dickson_construction

--
  Simen
January 02, 2014
On 01/01/2014 19:55, Joseph Rushton Wakeling wrote:
<snip>
> There are binary operations on complex numbers where the only sensible
> outcome seems to be non-integral real and imaginary parts. Addition,
> subtraction and multiplication are OK with integral types, but division
> really seems unpleasant to implement absent floating point,

Then why not just disable division if it's a non-float type, rather than preventing the whole complex template from being used with that type? This is like cutting off somebody's arm because they have a sore thumb.

Moreover, we have no way in the general case of determining whether T is an integral type, a library float-esque type, or (for example) a Galois field type.  So disabling it _just in case_ division doesn't work is crazy.  There must be better ways to do it.

> exponentiation even more so.

Exponentiation by a non-negative integer is straightforward.  So we should at least support this case for Gaussian integers.

<snip>
> However, I think relaxing the template constraints like this would best
> be done in the context of a library float-esque type (e.g. BigFloat)
> being implemented in Phobos, which could then be used to provide both
> proof-of-concept and the primary test case.

What do you mean by "in the context of", exactly?  Restricting it to some float-esque type that is in Phobos would still be overly restrictive.

>> Or even more exotically, use Complex!(Complex!real) to implement
>> hypercomplex numbers.
>
> Perhaps best implemented as a type in its own right?

It would, but removing the restriction would simplify the implementation of Hypercomplex(T) by enabling it to be a wrapper for Complex!(Complex!T).

Stewart.
January 02, 2014
On 01/01/2014 23:12, Lars T. Kyllingstad wrote:
<snip>
> /*  Makes Complex!(Complex!T) fold to Complex!T.
>
>      The rationale for this is that just like the real line is a
>      subspace of the complex plane, the complex plane is a subspace
>      of itself. Example of usage:

This doesn't seem to make sense - _any_ set is a subspace of itself.

But I suppose one way to look at it is that Complex!T means the minimal vector space containing x and i*x for all x in T.

>      ---
>      Complex!T addI(T)(T x)
>      {
>          return x + Complex!T(0.0, 1.0);
>      }
>      ---

So the point is to enable templated functions to accept a real or complex argument and return a complex number.

>      The above will work if T is both real and complex.
> */

T is "both real and complex"?  How is this possible?

Stewart.
January 02, 2014
On Thursday, 2 January 2014 at 11:37:22 UTC, Lars T. Kyllingstad wrote:
> I don't think we should worry too much about standards compliance. A library Complex type is quite different from a hardware floating-point type.

Are you sure?

Sometimes you need to translate an algorithm, you don't understand the inner workings of, from a codebase/cookbook. If std.complex differs from the most used c++/fortran implementations people will be confused, and you also end up having (machine translated) algorithm libraries each supplying their own complex type. Use 3 different libraries and you have to deal with 3 different complex types.

Floating point is rather sensitive to reordering of instructions, so I'd say you'll be better off mirroring one of the major existing implementations, otherwise accumulated discrepancies will be blamed on the language... A new tool that produce the same results as the old proven dinosaur tool look trustworthy. It makes you think that the conversion of your algorithms to the new tool was a success.
January 02, 2014
On 02/01/2014 15:38, Simen Kjærås wrote:
> On 2014-01-01 12:29, Stewart Gordon wrote:
>> Or even more exotically, use Complex!(Complex!real) to implement
>> hypercomplex numbers.
>
> For a more generic solution to this, see Cayley-Dickson construction[1]
> and my implementation of such:

Hypercomplex numbers are outside the scope of Cayley-Dickson construction.  The latter defines quaternions, octonions and so on.

Though it is confusing as there are at least two definitions of "hypercomplex number":
(a) a element of any number system that extends the complex numbers
(b) an element of a particular such number system, which is what I was using it to mean.

The hypercomplex numbers to which I was referring are as Fractint uses the term, and briefly described at
http://mathworld.wolfram.com/HypercomplexNumber.html

ij = k, i^2 = j^2 = -1, k = 1.  Multiplication is commutative.  But the distinctive thing about these is that they can be represented as an ordered pair of complex numbers, and then the complex multiplication formula works on them, hence my suggestion.

> https://github.com/Biotronic/Collectanea/blob/master/biotronic/CayleyDickson2.d
>
>
> It's nowhere near finished, has some bugs, and I haven't worked on it
> for at least a year, but it's an interesting way to create complex,
> hypercomplex, dual, and split-complex numbers (and combinations thereof).
<snip>

So it's a generalisation of Cayley-Dickson construction to produce a range of 2^n-dimensional algebras over the reals.  I'll have to look at it in more detail when I've more time.

Stewart.
January 02, 2014
On Thursday, 2 January 2014 at 18:12:56 UTC, Stewart Gordon wrote:
> Then why not just disable division if it's a non-float type, rather than preventing the whole complex template from being used with that type? This is like cutting off somebody's arm because they have a sore thumb.

Because that also seems to me to be an unpleasant option. A complex number implementation that fails to support ordinary arithmetic operations in all circumstances is pretty non-intuitive and will probably lead to unpleasant bugs in users' code.

> Moreover, we have no way in the general case of determining whether T is an integral type, a library float-esque type, or (for example) a Galois field type.  So disabling it _just in case_ division doesn't work is crazy.  There must be better ways to do it.

I don't follow your point here. We can constrain T however we see fit. The point isn't to have some perfect representation of every mathematical possibility, it's to have useful code that serves a good range of use-cases while being robust and easy to maintain. Restricting T to floating point types is a useful simplification here that has few costs in terms of expected use-cases.

> Exponentiation by a non-negative integer is straightforward.  So we should at least support this case for Gaussian integers.

Again, I don't see it as being useful to have a Complex!T whose support for binary operations may vary depending on the type of T.

> What do you mean by "in the context of", exactly?  Restricting it to some float-esque type that is in Phobos would still be overly restrictive.

No, I mean that until you have at least one actual float-esque type to test with, it is probably unwise to relax the template constraints that currently mandate built-in FP types.

> It would, but removing the restriction would simplify the implementation of Hypercomplex(T) by enabling it to be a wrapper for Complex!(Complex!T).

And complicate the implementation of Complex itself, for the sake of a likely very marginal special interest that could be supported quite well by an independent type.
January 02, 2014
On Thursday, 2 January 2014 at 19:32:32 UTC, Joseph Rushton Wakeling wrote:
> No, I mean that until you have at least one actual float-esque type to test with, it is probably unwise to relax the template constraints that currently mandate built-in FP types.

Fixed-point support might not be important on x86, but it has been needed and used for complex numbers in the past on other CPUs.
January 02, 2014
On Thursday, 2 January 2014 at 19:50:46 UTC, Ola Fosheim Grøstad wrote:
> Fixed-point support might not be important on x86, but it has been needed and used for complex numbers in the past on other CPUs.

I agree that should be supported. I just think that template constraints should only be relaxed in the context of concrete test cases, otherwise you're making promises you don't know you can keep.
January 02, 2014
On Thursday, 2 January 2014 at 19:58:27 UTC, Joseph Rushton Wakeling wrote:
> I agree that should be supported. I just think that template constraints should only be relaxed in the context of concrete test cases, otherwise you're making promises you don't know you can keep.

Yes, that is very true. Fixed-point is very sensitive to precision-loss so it requires some careful thinking, so it is definitively not a drop-in replacement type.

(It is actually useful on x86 too when you want reproducible results across binaries and hardware, such as in peer-to-peer syncing of computed state.)