Mod for complex arguments

[asmeurer]

How should Mod behave fo complex arguments?

[Shekharrajak]

Some links that may help:

1. http://math.stackexchange.com/questions/274694/modulo-complex-number
2. http://www.freemathhelp.com/forum/threads/76383-Modulo-of-complex-numbers
3. https://en.wikipedia.org/wiki/Gaussian_integer
4. http://fermatslasttheorem.blogspot.in/2005/06/division-algorithm-for-gaussian.html

[Shekharrajak]

Mod for complex number will help in these types of cases in solveset.

In [1]: solveset(exp(x) - sin(3),x)
Out[1]: {2⋅n⋅ⅈ⋅π + Mod(log(sin(3)), 2⋅ⅈ⋅π) | n ∊ ℤ}

{2⋅n⋅ⅈ⋅π + log(sin(3)) | n ∊ ℤ} ans looks better.

[siefkenj]

I think Mod should take an optional third argument integer_ring specifying what ring to take the modulus with respect to. Then mod(x, y, integer_ring=X) would return the representative of x in x.ambient_ring/(y*X). If x is a complex number, I think the default choice of ambient_ring should be the Gaussian integers ({a+b*i : a,b are integers})

[rpisarev]

I agree with @siefkenj. The sense of Mod between 2 complex numbers exists only in special rings, e.g. Gaussian integers or Eisenstein integers or even Kummer ring. We should have a norm for this numbers also. In this case we can go further: Euler's function phi(z), gcd(z_1, z_2) etc.

[asmeurer]

@Shekharrajak I don't think Mod is needed in that case. log is already assumed to be single-valued.

[Shekharrajak]

@asmeurer , It is coming from eval_imageset (canonical shift). I am not sure about this. (I think it is for like : expr = n*pi + 3*pi then n*pi + Mod(3*pi, pi) = > n*pi + 0 => n*pi . here match[a] = pi and match[b] = 3*pi )

[Shekharrajak]

I found a link.

Mod[m, n] == m - n * Floor[m/n]

where m, n is real or complex numbers.

for any complex number z

Floor(z) = Floor(Re(z)) + I*Floor(Im(z))

For example : Mod(-log(3), 2*I*pi ) = -log(3) - (2*I*pi) * Floor(-log(3)/(2*pi*I))

since Floor(-log(3)/(2*pi*I)) = 0. so -log(3) - (2*I*pi) * Floor(-log(3)/(2*pi*I)) = -log(3) .

Mod(-log(3), 2*I*pi ) = -log(3) is ans.
Is it correct ?

[asmeurer]

We do have that same definition for floor.

[asmeurer]

By that definition, larger real numbers mod 2*pi*I would be complex.

[Shekharrajak]

@asmeurer , For large number e.g. Mod(-log(3)*100000000, 2*pi*I)
Ans will be

(Mod[m, n] == m - n * Floor[m/n])

In [ ]: (-log(3)*100000000) - (2*pi*I)*floor(-log(3)*100000000/(2*pi*I))
Out[ ]: -100000000⋅log(3) + 34969914⋅π

same as wolframalpha

[Shekharrajak]

This method fails, when both m and n are complex numbers :
e.g.

Mod((2*I), ( 1+3*I))

In [23]: (2*I) - (1+3*I) * floor((2*I)/(1+3*I))
Out[23]: 2⋅ⅈ

but ans is -1 - I .

[jksuom]

If Mod should be defined for complex arguments, one should probably begin by considering its range of values for arbitrary complex arguments (and fixed modulus arguments).

On the real line, the range of Mod is an interval. It is typically a fundamental domain of a periodic function. The values of such a function in the fundamental domain determine its values elsewhere. The value at an arbitrary point is the same as the value at the Mod point in the fundamental domain.

In the complex plane there are many more possibilities for such fundamental domains. However, there are only two choices that are in general mathematical use.

• Periodic functions in the complex plane. For example, sin and exp. The fundamental domain is an infinite zone bounded by two parallel lines. The lines are usually taken orthogonal to the periods. The periods are the integral multiples of a fundamental period.
• Doubly periodic functions. The fundamental domain is a quadrilateral spanned by two fundamental periods that are linearly independent over the real numbers. The periods are the integral linear combinations of the fundamental periods.

[gschintgen]

I think the difficulty stems from the fact that there are (at least) two distinct notions at play that happen to be compatible in simple cases only:

1. When I read "modulo" I'm thinking about algebra, Euclidean division in Euclidean domains, etc. In the simplest case that would be computation in Z/nZ. But it could also mean the remainder of a Euclidean division over the Gaussian integers, as already indicated by @siefkenj and @rpisarev.

2. In an analytic context, there's the idea of "modulo" being used to reduce the argument of a periodic function to its fundamental domain. In simple cases (e.g. period = 2*pi), this can be considered a natural extension of the operation on Z/nZ. But just as well, we could be talking about doubly-periodic elliptic functions... Let's call this operation period_reduce for now.

Consider this example: (-3+9*I) % (1+I)

In the ring of Gaussian integers this would mean finding q, r (in C) such that -3+9*I = (1+I)*q + r where norm(r) <= norm(1+I), with norm(a+b*I) = a**2 + b**2 (no sqrt). In other words the remainder must lie in a disk centered at the origin. (E.g. q = 3+6*I, r = 0). If I'm not mistaken, the Wolfram definition (pointed to by @Shekharrajak) has the advantage that it returns a result that is coherent with Euclidean division of the Gaussian integers. As such it's a more or less a natural extension to C of this type of division, just as modulo a real number is an extension of the Euclidean division over Z.

But what if the user actually intended to period_reduce an argument to the (contrived/meaningless?) complex function f = sin(2*pi*x/(1+I)) which is periodic with a period of 1+I? In that case, the answer would be -3+9*I - n*(1+I) where n is a suitable integer (here: n = 3) so that the result lies in-between the two parallel lines of the fundamental domain, as explained by @jksuom.

Given all of the above, I'd say there are two possible choices:

1. Out of precaution Mod could raise an exception (ValueError?) if called with any non-real arguments. That would force the user to explicitly use either some algebraic method or a (yet to be implemented) period_reduce.
2. Use the "Wolfram definition", but only if it's determined to be compatible with Sympy's current Mod which itself is defined to be compatible with Python's conventions...

Anyway, in the above-mentioned (unnecessary) Mod(log(...), 2*pi*I) I'd expect the library code to either extract the imaginary part and reduce it "modulo 2*pi" or to call some period_reduce(..., 2*pi*I) method. In my opinion relying on modulo to period-reduce a complex number is simply wrong. See #11391 (comment)

[gschintgen]

Just a short observation using 1.5:

In [1]: Mod(2+3*I, 5-4*I)
Out[1]: (2 + 3⋅ⅈ) mod 5 - 4⋅ⅈ

In [2]: Mod(2+3*I, 5-4*I).rewrite(floor)
Out[2]: 7 - ⅈ

Of particular note is that

• Wolfram gives -2-2i.
• If any alternative definition to Sympy's current a-b*floor(a/b) (used by rewrite) were to be retained based upon the previous discussion in this thread, there'd be an incompatibility between plain Mod(complex, complex) and Mod(complex, complex).rewrite(floor).

[gschintgen]

Hmm, Mod does in fact act on (some) non-reals:

In [1]: (3+5*I) % 2
Out[1]: (1 + 5⋅ⅈ) mod 2

This auto-simplification is unexpected, especially with the mod 2 still present. It's also incompatible with the formula used by .rewrite(floor).

[gschintgen]

More strangeness, this time % acts on the imaginary part, contrary to the previous example:

In [4]: (3*I) % 2
Out[4]: ⅈ mod 2

[asmeurer]

The floor definition seems like a good fundamental thing to keep. That means that a == a//b*b + a%b (in SymPy, a//b is floor(a/b)). floor already generalizes to complex numbers. If I understood the comments above, Wolfram also uses this definition, assuming they use the same floor definition as SymPy.

This definition is used to match the Python definition of mod for negative arguments (which doesn't match other languages, see https://stackoverflow.com/questions/3883004/the-modulo-operation-on-negative-numbers-in-python).

[gschintgen]

The floor definition seems like a good fundamental thing to keep.

This is going to be a long post...

I actually wanted to just say "Agreed, let's keep the floor definition, just like Wolfram", but then something came up while playing around with WolframAlpha. Here are my findings:

• Sympy-floor is specified to take "the floor of the real and imaginary parts separately". This corresponds to the Wolfram language: "Floor applies separately to real and imaginary parts of complex numbers." I did a few trial computations and did not find any obvious difference.

• Wolfram's Mod is defined as m-n*floor(m/n) (according to functions.wolfram.com) or m-n*Quotient[m, n] (according to reference.wolfram.com). Here Quotient refers to floor(m/n) (according to functions.wolfram.com) but the formulation seems to restrict this to integers. And this weaseling around is there for a reason:

Quotient[-3+i, 2-i] == -1
Floor[(-3+i)/(2-i)] == -2-i
SymPy: (-3+I)//(2-I) == -2-I == floor((-3+I)/(2-I))

The numerical value of the quotient is -1.4-0.2*I, so it's as if Quotient were defined to be round(re(a/b))+I*round(ima/b)). This is actually not that bad, since rounding instead of floor leads to a correct way to compute a/the remainder of a Euclidean division over the Gaussian integers. (See Wikipedia)

Obviously this is not compatible with the standard floor-based definition on the real line.
But it's even stranger:

Quotient[-5, 2]        == -3        (floor-based definition)
Quotient[5, 2]         == 2         (floor-based definition)
Quotient[-7, 2]        == -4        (floor-based definition)
Quotient[7, 2]         == 3         (floor-based definition)
Quotient[-5+0.01*i, 2] == -2        (round-to-even)
Quotient[5+0.01*i, 2]  == 2         (round-to-even)
Quotient[-7+0.01*i, 2] == -4        (round-to-even)
Quotient[7+0.01*i, 2]  == 4         (round-to-even)

And this muddled definition/implementation then leads to the following bug:
For Quotient[-3+i, 2-i] wolframalpha lists as "Alternative representation" the floor-based definition, including the copy&pastable string Quotient[-3 + I, 2 - I] == Floor[(-3 + I)/(2 - I)] which, when given as input, returns False!

My best guess is that Wolfram did recognize that the floor-based definition is not really mathematical useful on the complex plane, so they decided to concoct some Frankensteinian creation that

• corresponds to the floor-definition over the reals,
• gives a correct quotient/remainder over the Gaussian integers as a Euclidean domain,
• uses round-to-even for unknown reasons,
• respects n*Quotient[m,n]+Mod[m,n] always being equal to m,
• leads to bugs due to its peculiar definition and sub-par documentation...

To sum it up:

• SymPy is already incompatible with Wolfram when it comes to integer division: // vs. Quotient.
• Mod over Gaussian integers is incompatible with floor-based definition, do not try to reconcile them.

IMHO there are two options:
A) Since the floor-based definition on the complex plane is not that useful mathematically (or am I missing something?), raise an exception if called with non-real arguments and document this design decision.
B) Use the floor-based definition, just so that the Mod function is technically defined over the whole complex plane and document it.

In any case, especially B), Mod's docstring should mention the caveats:

• do not use Mod a.k.a % to reduce modulo some number if any of the two arguments may be non-real, including counter-examples (e.g. (1+3*pi*I) % (2*pi*I) is not 1+pi*I; Mod(5*pi + 10*I, 2*pi) is not pi + 10*I)
• do not mistake Mod for a Euclidean division remainder over the Gaussian integers. Counter-example: Mod(1-4*I, 3+I) == 2+3*I but norm(2+3*I)>norm(3+I).

[asmeurer]

I agree whatever we need to do we definitely need to document it. The documentation should also include the proper ways to do both of those things. We sometimes fall short in documenting reasons on why we define certain functions in certain ways, especially regarding extensions of functions normally defined on integers to complex numbers. It's still on my todo list to document the discussion from #14577.

do not use Mod a.k.a % to reduce modulo some number if any of the two arguments may be non-real, including counter-examples (e.g. (1+3piI) % (2piI) is not 1+piI; Mod(5pi + 10I, 2pi) is not pi + 10*I)

Did we discuss above a definition of taking mod of the real and imaginary parts separately? That would do that. There's also a - re(a//b)*re(b) - im(a//b)*im(b)*I which works for the above examples although I'm not sure how to interpret it when both numbers have nonzero real and imaginary parts (it likely isn't meaningful).

I personally wouldn't be surprised that Mod doesn't work over Z[i], but I would expect it to be useful for reducing complex modulus.

[gschintgen]

Did we discuss above a definition of taking mod of the real and imaginary parts separately?

No, but I can see its appeal: It's simple and pragmatic enough and could be documented along the lines "Mod is intended for reals only but for convenience the following behavior for complex arguments has been implemented (...) The rewrite(floor) method is only equivalent over the reals."

That would do that. There's also a - re(a//b)*re(b) - im(a//b)*im(b)*I

for a=10+7*I and b=2*I that would give 10+17*I.

[asmeurer]

We can change rewrite(floor) if we go with that. It still reduces to the formula for both arguments real.

[gschintgen]

We can change rewrite(floor) if we go with that. It still reduces to the formula for both arguments real.

Oh, that's right of course. Even better.