Syntax for complex numbers

[Atomotron]

Idris 2 supports complex numbers, but lacks an easy way to write them. Presently, the only way I know to get the imaginary unit i is with,
sqrt (-1)
or,
log (-1) / pi
Since stringification will happily print 0.0 + 1.0i after evaluating log (-1) / pi, and since the Double (cast "0.0+1.0i") will return successfully, it stands to reason that it should work the same way in source code. Most other types have a bijection between parsing and printing, and it stands to reason that the complex part of the double type should also.

Of course it is possible to write expressions that evaluate to complex numbers, but that lacks the nice property of the printed form of types being valid source code. Idris would be best if it only had one Double parsing behavior, and not two subtly distinct behaviors, one for parsing source, and the other for cast.

[gallais]

Wait, does it?

The type of sqrt is Double -> Double so I don't see why it
should return a complex number (it's not even well-defined on -1).

This seems more of a bug due to a weird thing the backend decides to
do with Double?

[Atomotron]

Here are a few lines from a recent REPL session on the latest master branch:

Main> :let x : Double
Main> :let x = sqrt(-1)
Main> x
0.0+1.0i

I think it is weird that the Double type can do this, especially because of its correspondence with doubles in the C API. However I love that Idris supports complex numbers natively, because I write a lot of programs that use them. :)

(it's not even well-defined on -1).

sqrt is well-defined on negative numbers within the algebra of complex numbers. A branch cut must be chosen, but since the same branch cut is chosen almost all the time, that's not a very important technicality. Idris supports algebra over complex numbers quite well, for example:

Main> sqrt(sqrt(-1))
0.7071067811865476+0.7071067811865475i

All of the other basic math functions work over the complex numbers as well, for example log and exp.

[gallais]

I think this is an undefined behaviour. The expected behaviour is probably
that it should raise an exception if called with a negative number.

Chez scheme seems to indeed specify that it'll try to make the best of it
but I wouldn't be surprised if you can break things by using a different
backend: https://www.scheme.com/tspl4/objects.html#./objects:s127

[gallais]

There we go:

idris2 --cg node --no-banner
Main> :let main : IO ()
Main> :let main = print (sqrt (-1))
Main> :exec main
NaNMain> 

vs

idris2 --no-banner
Main> :let main : IO ()
Main> :let main = print (sqrt (-1))
Main> :exec main
0.0+1.0iMain> 

[Atomotron]

Do you think it would be better to return NaN like 0.0/0.0 does, or to violate totality? I think it might be better to keep totality, if for nothing else than to ease totality proofs. (It would be difficult to prove that a double would be in the range of acos, given that it's delimited by transcendental numbers. I am not really sure how proofs about doubles would work in Idris, anyway. I have only seen proofs about Nats.)

If acos, log, sqrt and the other basic math functions weren't total, I would foresee either partiality cascading through programs, or a bunch of "just trust me" compiler directives.

[stefan-hoeck]

Returning NaN in case of 0.0/0.0 does not violate totality. NaN is a valid Double according to the IEEE standard (though these constants are right now not officially supported by Idris2, see #29) and the program does not crash when returning a NaN result. It should be possible to refine Double in your own code (by using Data.DPair.Subset for instance) in such a way that it does not include NaN or Infinity.

Concerning proofs about floating point numbers: Alas, these will be very hard to do properly. For instance, even the following most basic axiom does not hold for Double (even without the above mentioned constants): x,y > 0 => x * y > 0.

Edit: Neither does the law of distributivity:

Main> 10 * (1.123 + 1.0e-16)
11.23
Main> 10 * 1.123 + 10 * 1.0e-16
11.230000000000002

[Atomotron]

Returning NaN in case of 0.0/0.0 does not violate totality.

Yes, I agree. I think it would be better to return NaN for expressions like sqrt(-1) than it would be to crash on a misapplied partial function, both of which are better than undefined behavior. I don't think Maybe Double would be a good return type, because all of those checks could really slow down numeric code, even if the programmer would prefer to let IEEE propagation of NaNs handle it. Like you said, a programmer intent on having the type system check their numeric code could start by defining their own safe double type, although you correctly observe that proving anything will be difficult.

[ohad]

I support making sqrt partial, and adding a safe version in the library for use in proofs (using proof reflection):

import Syntax.WithProof

sqrtSafe : (d : Double) -> {auto 0 prf : d >= 0 = True} -> Double
sqrtSafe d = assert_total (sqrt d)

One : Double
One = sqrtSafe 1

{-
Error : Double
Error = sqrtSafe (-1)  
-- While processing right hand side of Error. Can't find an implementation for False = True.
-}

Complicated : Double -> Double
Complicated x = case @@(x >= 0) of
  (True ** prf) => sqrtSafe x
  (False** _  ) => 0
  
Zero : Double
Zero = Complicated (-1)

[Atomotron]

sqrt on doubles should be partial iff division on doubles is partial. As far as I can tell the same reasons apply to both, (Fidelity to IEEE standards, similarity with other languages, fidelity to typical definitions...), so whatever the right answer is, it would be the same answer in both cases. It also seems like it would be bad to have some operations on doubles return NaN for mathematically undefined cases while others crash.

[ShinKage]

IMHO as @Atomotron is saying, consistency is the most important factor. Double should be treated as a primitive implementing the IEEE standard with its sentinel values such that all operation are total even if mathematically undefined. Library abstractions should, instead, be in charge of providing a more safe abstraction on top of Double which may be closer to the mathematical approach. Now that we also have double literal overloading, it would be also an ergonomic approach.

[Atomotron]

Another point to consider is that if there will ever be a codegen targeted at high performance, mirroring IEEE (i.e. whatever Intel believes should happen) will be the only reasonable option - as putting checks after every sqrtsd to conform to a language specification that demands crashing would not be appropriate.

My understanding is not perfect, but I think that although hardware FP exception handling can be configured, it involves things like getting "SIGFPE" signals from the operating system[0]. We probably do not want to be getting signals from the operating system, and if we don't want to be emitting branches after every sqrt instruction, that leaves NaN as the remaining option.

Lua, Javascript, C, and probably LLVM IR (since x86 does it) backends will all see simpler output if NaN is returned for out-of-domain calls to transcendental functions.

Edit: I found this (flsqrt fl) alternate sqrt function in the Chez scheme documentation that always returns flonums. Maybe fixing this issue would be as simple as having the default codegen use these.

[0] https://www.jayconrod.com/posts/33/trapping-floating-point-exceptions-in-linux

[gallais]

Yes, I think you're right: consistent, well-specified behaviour for the builtins
plus a (choice of) library layers offering a safer interface (something à la Flocq
would be a natural candidate).