Extend floor (round, etc) to Complex numbers

[danilo-bc]

I was recently porting a Matlab example to julia and I realized there's no direct equivalence.

I extended the floor function in my code snippet by doing:

import Base.floor
floor(x::Complex) = floor(real(x)::Real) + 1im*floor(imag(x)::Real)

I'm not sure if this would be a welcome addition, but it would surely make Julia code more similar to other well-used scientific software.

[oscardssmith]

If there is a concensus among other languages that this is the right behavior, it seems reasonable to me. If so, it would be a good first PR.

[Seelengrab]

Wolfram Alpha has this definition as well, but it doesn't seem that clear cut, as floor requires some interpretation of an inequality, which we don't have for complex numbers (and neither does Wolfram Alpha).

Existing implementations seem to project onto gaussian integers (like in the OP here), but there seem to be at least some issues with that in regards to some identities. Python doesn't have a complex floor, as far as I can tell (does numpy?). Not sure what the right thing to do is 🤷‍♂️

[Seelengrab]

There's also a list of interesting requirements for a floor(x::Complex) function based on domain extension done by Eugene McDonnel for the J language (and similarly done for Dyalog APL and NARS2000), reproduced (and translated to julia) here for convenience:

1. Existence: Every number has a floor.
2. Uniqueness: Every number has only one floor.
3. Fractionality: The magnitude of the difference of a number and its floor shall be less than one. This property must be satisfied to guarantee that remainders are less in magnitude than divisors. It may be called the fundamental property of the floor function.
4. Integrity: The floor of a number is an integer. (note: gaussian integer is meant here, i.e. a complex number with both and real and imaginary part being an integer from Z)
5. Convexity: If g is the floor of the numbers z and w , then it is also the floor of all numbers on the line segment between z and w .
6. Integer Translation: For c a complex integer, c+floor(z) == floor(c+z).
7. Compatability: The complex floor function is compatible with the real floor function. Furthermore, its action on purely imaginary numbers is similar to the action of the real floor function on real numbers. In particular, real(floor(z)) ≤ real(ceil(z)) and imag(floor(z)) ≤ imag(ceil(z)) .

The naive implementation of just flooring components violates Fractionality, as noted in this codegolf (which also has a version with allegedly minimal rounding required).

It should be noted that if we go the APL route, we can similarly define ceil(::Complex) and apparently get gcd(::Complex) for free.

[danilo-bc]

Great insights! In this case, I feel that it's better not to have a floor(x::Complex) at all.

This is because people may port their codes from somewhere else, have it not raise an error, and not find functional bugs.

For my application result = floor(real(x)) + 1im*floor(imag(x)) seems to be the proper solution, instead of having the extended version.

[Seelengrab]

I don't know if the APL version is "more correct" - the identity floor(sqrt(floor(z))) == floor(sqrt(z)) also doesn't hold for the APL version (due to us choosing one root over another, I believe?). Kind of a bummer, and I don't own "Concrete Mathematics: A Foundation for Computer Science" by Graham, Knuth and Patashnik so I'm afraid I can't investigate further for maybe a better version that would allow for that and more identities to hold.

[Seelengrab]

If others want to play around with either version, here are both together (MIT, free to use, yada yada):

naiveFloor(z::Complex) = floor(real(z)) + floor(imag(z))im

function aplFloor(z::Complex{T}) where T
   r = real(z)
   i = imag(z)
   b = floor(r) + floor(i)im
   # not sure whether this should be `rem` or `mod`..
   x = mod(r, one(T)) 
   y = mod(i, one(T))
   if one(T) > x+y
      return b
   else
      return b + (T(x ≥ y) + T(x < y)im)
   end
end

For my application result = floor(real(x)) + 1im*floor(imag(x)) seems to be the proper solution, instead of having the extended version.

If it works for reproducibility in your port, that's fine! Best to be aware of the caveats though.

[eschnett]

@Seelengrab Does the "magnitude" in the fractionality constraint need to be the "usual" magnitude, or could it be e.g. based on the infinity-norm? In this case, magnitude(x::Complex) = max(abs(real(x)), abs(imag(x))), and the fractionality constraint holds. The infinity-norm is arguably more suited for Gaussian integers, which are compatible with the result type Complex{Int}. (What one loses this way is the radius/phase representation of complex numbers.)

This definition extends nicely to other vector spaces. If one treats complex numbers as a vector space over reals, this is a natural extension.

[Seelengrab]

Does the "magnitude" in the fractionality constraint need to be the "usual" magnitude, or could it be e.g. based on the infinity-norm?

All references to magnitude in the sources I cited refer back to Gauss in the end, so it's the "usual" meaning, not any special norm. The difficulty is exactly that you can't really just pick a norm you like and say "this is the general norm we use", because it may not be suited for all purposes (as encountered in this thread). See the linked APL paper for the references and background in question - there's also some nice graphs.

In this case, magnitude(x::Complex) = max(abs(real(x)), abs(imag(x))), and the fractionality constraint holds.

Keep in mind that the magnitude itself is irrelevant for floor (aside from requiring the magnitude of the floor to be the smallest of all the numbers flooring to that) - choosing a particular definition for magnitude effectively means nudging the definitions around to make that fit. It's also a bit problematic to choose this, since norm(::Complex) exists already and is not that, so choosing that definition for this case only would introduce an inconsistency.

---

I don't think choosing that definition for magnitude is necessarily better or worse, it's just different. I think that one consequence of requiring that would be that all numbers in the unit square map to the origin, tiling along integer multiples in each dimension, which I think just comes back to the naive flooring of components, no?

[adienes]

A priori I don't have a strong opinion over whether complex(floor.(reim(z))...) is an appropriate floor function for complex numbers.

However, due to

julia> isinteger(1 + 1im)
false

I don't think it should be added now, unless breaking changes to isinteger are made to return true for Gaussian ints

[LilithHafner]

In Julia 1.0 we have

julia> round(1.1 + 2.9im, RoundDown)
1.0 + 2.0im

help?> round
...
  Return the nearest integral value of the same type as the complex-valued `z` to `z`...

[Seelengrab]

That method can take two rounding modes though, not just one:

[...], breaking ties using the specified RoundingModes. The first RoundingMode is used for rounding
the real components while the second is used for rounding the imaginary components.

thereby sidestepping the issue of consistency. This issue is about floor(::Complex) by itself not having a good default (and arguably, the single argument round(::Complex) also shouldn't exist, for the same reasons).

[Seelengrab]

Reopening since it was not clear cut why this should be implemented as a fall back.