January 02, 2014
On 02/01/2014 19:32, Joseph Rushton Wakeling wrote:
> 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.

The compiler rejecting the code is the most pleasant bug that's possible IMO.

>> 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.

Not being overly restrictive in what types you will allow it to be used with is an important part of serving a good range of use cases.

<snip>
>> 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.

I'm sure we have a small handful already.  We just need to find them. For instance, given the time I could probably dig up my rational number implementation and update it to current D.

>> 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.

How would it complicate the implementation?  Removing the undocumented rule whereby Complex!(Complex!T) folds to Complex!T would be a slight simplification.

Maybe the right course of action is to have a parameter to the template that suppresses the type restriction and the folding rule, so that people who want to use it on a type that might not work properly can do so.  This would be a relatively small change.

Stewart.
January 02, 2014
On Sunday, 24 November 2013 at 17:35:34 UTC, Shammah Chancellor wrote:
> On 2013-11-24 15:50:46 +0000, Joseph Rushton Wakeling said:
>
>> On Saturday, 23 November 2013 at 15:13:22 UTC, Shammah Chancellor wrote:
>>> I disagree.  I was using them for physics simulations.   They are very useful for the computational physics community.   Just because most people are still using FORTRAN does not mean they won't switch eventually.
>> 
>> Would it cause you any particular disadvantage to use the library std.complex rather than the built-in complex type?
>
> It would if the they don't work correctly.  There needs to be an Imaginary type and some proper operations between complex and imaginary types.  That doesn't seem to be the case currently.  I personally think having the built-in type is very helpful.  However, I can understand from a language perspective that having "i" around is hard for the parser.
>
> Also, the argument "If complex/imaginary is built-in, why not have quaterions also" seems to imply that it should be a library type.
>
> -Shammah

You can have im!5.0 in a library type, to me it sounds good enough.
January 02, 2014
On Thursday, 2 January 2014 at 20:23:40 UTC, Stewart Gordon wrote:
> The compiler rejecting the code is the most pleasant bug that's possible IMO.

It's still unpleasant relative to the intuitive expectation that numerical types should support all basic numerical operations.

> Not being overly restrictive in what types you will allow it to be used with is an important part of serving a good range of use cases.

It's a question of balance: breadth of support vs. ease of implementation and maintenance.  Some use cases may be better served in other ways than rolling them into one "catch-all" type.

> I'm sure we have a small handful already.  We just need to find them. For instance, given the time I could probably dig up my rational number implementation and update it to current D.

I'd really like to see that anyway, for its own sake; I've been working on a std.rational based on David Simcha's work, it would be good to compare.

> How would it complicate the implementation?  Removing the undocumented rule whereby Complex!(Complex!T) folds to Complex!T would be a slight simplification.

It would involve non-trivial revisions to the internals of Complex.

If you want to submit a pull request supporting this I'm sure it'll be considered carefully, but it is not something to be done casually.

> Maybe the right course of action is to have a parameter to the template that suppresses the type restriction and the folding rule, so that people who want to use it on a type that might not work properly can do so.  This would be a relatively small change.

I do not know how familiar you are with the internals of std.complex.Complex but I think it would be a good idea to go through them in some detail before asserting that changes like this would be "relatively small".

January 02, 2014
On Thursday, 2 January 2014 at 01:17:54 UTC, Joseph Rushton Wakeling wrote:
>    * There are calculations involving "pair of reals" implementations of complex
>      numbers that generate spuriously incorrect results. They're corner cases,
>      but they exist.
>
>    * A purely imaginary type allows those calculations to be performed
>      correctly. It also allows more precise calculations in general for cases
>      where the real part of a complex number is known to be 0.

So we have three choices:

1. We add a dedicated Imaginary type, and leave Complex more or less the way it is. The pros are that it provides a way to represent imaginary numbers "correctly", while keeping Complex simple and performant. The cons are that Complex still gives wrong (or at least somewhat unexpected) results in the aforementioned corner cases, and that we add a new type with extremely limited usability.

2. Special-case Complex for imaginary numbers. The pros are that it solves the problems Imaginary was intended to solve, and we don't need a new type. The cons are that the Complex implementation becomes more complex (haha) and less performant, and that we actually change the semantics of IEEE floating-point "zero" somewhat.

3. Leave std.complex as-is, and make sure people know what the problematic cases are.

They all suck. I don't know what is the lesser of three evils here, but I too am starting to lean towards 1.  I'm probably going to continue playing devil's advocate, though. ;)


>    * You can do calculations involving purely real-valued numbers and complex
>      numbers and not run into the same issues, because purely real values are
>      supported. So you should be able to do the same with purely imaginary
>      numbers.

That argument is fallacious. Imaginary numbers are quite different from real numbers.


> I should add that as far as I'm concerned what I want is simply to find the best possible way to represent and use purely imaginary numbers in D. I don't mind if, having done that, the code winds up being rejected, as long as the exercise proves useful.

I think it's great that you're doing this. :)
January 02, 2014
On 02/01/14 22:08, QAston wrote:
> You can have im!5.0 in a library type, to me it sounds good enough.

As currently written, it'll be imaginary(5.0) or Imaginary!double(5.0) -- sorry for the verbosity, but "im" seems to me too likely to wind up being used as a variable name ... ;-)

January 02, 2014
On 02/01/14 23:26, Lars T. Kyllingstad wrote:
> 2. Special-case Complex for imaginary numbers. The pros are that it solves the
> problems Imaginary was intended to solve, and we don't need a new type. The cons
> are that the Complex implementation becomes more complex (haha) and less
> performant, and that we actually change the semantics of IEEE floating-point
> "zero" somewhat.

Well, what I was mulling over was a Complex type that consists of two FP values (re and im) and two bools (let's call them hasRe, hasIm).

The point of this is that if you assign to Complex from a regular real-valued number, then the value of re gets set, im gets set to zero, hasRe to true, and hasIm to false.  So, you have a flag inside that says, "Hey, only worry about your real part."

Then you could define an imaginary() helper function that creates a Complex number with re = 0.0, im set to whatever it needs to be, hasRe to false, and hasIm to true.  So again, you can benefit from a Complex that _knows_ it's purely imaginary.

Moreover, it should be possible to preserve that kind of knowledge under certain operations.  e.g. if you multiply two Complex numbers which both have hasRe == false, you can _know_ that you're going to wind up with a number where hasRe == true and hasIm == false.

OTOH when you perform certain other operations, you might end up with a zero re or im, but also with hasRe or hasIm to true, which basically reflects your _certainty_ that the real or imaginary part is zero.

It might be interesting to try out, and it seems nicer to me than if ()'s based on whether re == 0 or im == 0, but on the other hand it enlarges the size of a Complex type and will result in a performance hit, and will break compatibility with C/C++ Complex types.

> They all suck. I don't know what is the lesser of three evils here, but I too am
> starting to lean towards 1.  I'm probably going to continue playing devil's
> advocate, though. ;)

Please do.  It's fun, and it also helps to improve things. :-)

>>    * You can do calculations involving purely real-valued numbers and complex
>>      numbers and not run into the same issues, because purely real values are
>>      supported. So you should be able to do the same with purely imaginary
>>      numbers.
>
> That argument is fallacious. Imaginary numbers are quite different from real
> numbers.

Can you expand on that?

> I think it's great that you're doing this. :)

Thank you! :-)

January 02, 2014
On 02/01/2014 21:10, Joseph Rushton Wakeling wrote:
> On Thursday, 2 January 2014 at 20:23:40 UTC, Stewart Gordon wrote:
<snip>
>> How would it complicate the implementation?  Removing the undocumented
>> rule whereby Complex!(Complex!T) folds to Complex!T would be a slight
>> simplification.
>
> It would involve non-trivial revisions to the internals of Complex.

Please be specific.

<snip>
>> Maybe the right course of action is to have a parameter to the
>> template that suppresses the type restriction and the folding rule, so
>> that people who want to use it on a type that might not work properly
>> can do so.  This would be a relatively small change.
>
> I do not know how familiar you are with the internals of
> std.complex.Complex but I think it would be a good idea to go through
> them in some detail before asserting that changes like this would be
> "relatively small".

Oh yes, some of the templated functions that take or return a complex would need this extra parameter adding.

Another way it could be done is to have a Complex template that implements these rules and which programmers would normally use, and have this aliasing ComplexImpl which actually provides the implementation and which programmers can use directly if they want to bypass the restrictions.  The difficulty is documenting it in a way that would make sense to normal users....

Stewart.
January 02, 2014
On Thursday, 2 January 2014 at 18:37:36 UTC, Ola Fosheim Grøstad wrote:
> 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?

Not at all. ;)

I just think we should keep in mind why FP semantics are defined the way they are. Take 0.0*inf, for example.  As has been mentioned, 0.0 may represent a positive real number arbitrarily close to zero, and inf may represent an arbitrarily large real number. The product of these is ill-defined, and hence represented by a NaN.

0.0+1.0i, on the other hand, represents a number which is arbitrarily close to i. Multiplying it with a very large real number gives you a number which has a very large imaginary part, but which is arbitrarily close to the imaginary axis, i.e. 0.0 + inf i. I think this is very consistent with FP semantics, and may be worth making a special case in std.complex.Complex.


> 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.

I agree, but there is also a lot to be said for not repeating old mistakes, if we deem them to be such.
January 03, 2014
On Thursday, 2 January 2014 at 23:26:48 UTC, Joseph Rushton Wakeling wrote:
> On 02/01/14 23:26, Lars T. Kyllingstad wrote:
>>>   * You can do calculations involving purely real-valued numbers and complex
>>>     numbers and not run into the same issues, because purely real values are
>>>     supported. So you should be able to do the same with purely imaginary
>>>     numbers.
>>
>> That argument is fallacious. Imaginary numbers are quite different from real
>> numbers.
>
> Can you expand on that?

Mathematically, the real numbers are a ring, whereas the purely imaginary numbers are not. (I've been out of the abstract algebra game for a couple of years now, so please arrest me if I've remembered the terminology wrongly.)  What it boils down to is that they are not closed under multiplication, which gives them radically different properties - or lack thereof.
January 03, 2014
On 03/01/14 00:33, Stewart Gordon wrote:
> Please be specific.

You know, I'm starting to find your tone irritating.  You are the one who's asking for functionality that goes beyond any Complex implementation that I'm aware of in any other language, and claiming that these things would be trivial to implement.  I would expect a person who claims with confidence that something is trivial, to actually know the internals of the code well enough to understand what parts of it would need to be modified.  On the basis of what you've written, I have no reason to believe that you do.

There are numerous places inside current std.complex.Complex where temporary values are used mid-calculation.  Those are all of type FPTemporary (which in practice means real).  So, to handle library types (whether library floating-point types such as a BigFloat implementation, or a Complex!T so as to support hypercomplex numbers) you'd either have to special-case those functions or you'd have to provide an alternative Temporary template to handle the temporary internal values in the different cases.

You'd also need to address questions of closure under operations (already an issue for the Imaginary type), and precision-related issues -- see e.g. the remarks by Craig Dillabaugh and Ola Fosheim Grøstad elsewhere in this thread.

While it could be done, I don't think it would be simple to get right and I would regard it as a major revision of std.complex.  I'd also be concerned that generalizing Complex in this way might lead to performance regressions for the standard Complex!T use-cases, due to the more complicated templates involved.

> Oh yes, some of the templated functions that take or return a complex would need
> this extra parameter adding.

No.  The internals would need significant rewriting, for reasons I've already elaborated.

> Another way it could be done is to have a Complex template that implements these
> rules and which programmers would normally use, and have this aliasing
> ComplexImpl which actually provides the implementation and which programmers can
> use directly if they want to bypass the restrictions.  The difficulty is
> documenting it in a way that would make sense to normal users....

If you want these things done in these ways, then I think the onus is on you to provide code which works and doesn't damage existing std.complex functionality.

I don't personally see any rationale for implementing hypercomplex numbers as a specialization of Complex given that they can be just as well implemented as a type in their own right.