What to do with mod?

[dpsanders]

How should we deal with functions like mod acting on intervals?

E.g.

• mod(0.8..1.2, 1) should give (0..0.2) ∪ (0.8..1);
  mod(0.8..2.2, 1) should give 0..1?

One solution is to define an interval union type; see e.g.
https://link.springer.com/article/10.1007/s10543-016-0632-y

Really the correct solution is probably to use contractors.

[lbenet]

Speaking about Interval, I think it should behave as tan. For the specific example you give, I think it should return 0..1, which corresponds to the hull.

A complementary possibility is to have mod(a::Interval) returning the hull, and also implement extended_mod(a::Interval), among other examples (e.g. extended_inv(a::Interval)).

extended_mod(a) would return a Vector{Interval) (or perhaps NTuple{N,Interval}) containing the list of (connected) intervals, as a way of representing the union in the sense you use above. If the interval a does not contain the discontinuity, the result would correspond to [mod(a)].

Perhaps, this is what you mean by using contractors.

[lbenet]

In #145 I have implemented something that is equivalent to mod(0.8..2.2, 1) --> 0..1.

We could also implement some extended_mod which would return a tuple of intervals (two or three, this we should discuss) which somehow mimics extended_div.

[dpsanders]

Great, thanks. Why would there be three intervals? Is this version of mod (which just returns 0..1) useful? I tried it with root finding and it didn't seem to be too useful.

…

On Thu, 3 May 2018, 3:01 pm Luis Benet, ***@***.***> wrote: In #145 <#145> I have implemented something that is equivalent to mod(0.8..2.2, 1) --> 0..1. We could also implement some extended_mod which would return a tuple of intervals (two or three, this we should discuss) which somehow mimics extended_div. — You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub <#129 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AALtTniH72rkhmw7slmA8eqx0WR-1EQJks5tu1QFgaJpZM4TXHeO> .

[lbenet]

Why would there be three intervals?

Considering your example, mod(0.8..2.2, 1), I would say that extended_mod would return (0.8..1.0) ∪ (0.0..1.0) ∪(0.0..0.2), if we consider each segment.

Is this version of mod (which just returns 0..1) useful?

This may be useful for things like find_quadrants; for root-finding methods, we would need extended_mod. I'll push a new commit dealing with extended_mod later.

[lbenet]

Just pushed another commit to #145 with extended_mod implemented. Can you retest your example with root finding, using now extended_mod? Maybe you could post it so I can play with it.

[dpsanders]

The problem is that all of the root finding functions assume that a single box in gives a single box out. Basically extended mod etc lead to interval unions, which is not yet implemented.

…

On Fri, 4 May 2018, 11:49 am Luis Benet, ***@***.***> wrote: Just pushed another commit to #145 <#145> with extended_mod implemented. Can you retest your example with root finding, using now extended_mod? Maybe you could post it so I can play with it. — You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub <#129 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AALtTnG9ekSMND3_12RR-cNS2dhdwRIDks5tvHihgaJpZM4TXHeO> .

[lbenet]

I'm not sure I understand your point.

I agree that extended_mod returns an interval union (as a tuple), which was precisely the suggestion made above. It mimics what Newton method does by extending the division, thus returning a tuple of disjoint intervals, if the derivative contains a zero. Similarly, extended_mod returns a tuple; each interval corresponds to the range of mod when the domain has no discontinuities.

[dpsanders]

I just mean that the current versions of the algorithms can't deal with interval unions.

[Datseris]

Is there any update here, or at least instructions of how to handle mods? I have a simple dynamical system that involves the mod operation

function dissipative_standard_map_rule(u, p, n)
    x, y = u
    (; ν, f₀) = p
    s = x + y
    xn = mod2pi(s)
    yn = (1 - ν)*y + f₀*(sin(s))
    return SVector(xn, yn)
end

and I would like to find its fixed points using interval arithmetic (root finding) via utilizing Newton's method (I am using ChaosTools.fixedpoints). However, inevitably, the mod is called on an interval and I get

ERROR: MethodError: no method matching rem2pi(::IntervalArithmetic.Interval{Float64}, ::RoundingMode{:Down})

Since I am anyways only searching for roots for x bounded in [0, 2pi], would you say it is correct to simply remove the mod operation from the above function?

[lbenet]

I guess the work you are looking for is in #145, but somehow we got distracted with other stuff...

[Datseris]

Seems that #525 is much more recent and perhaps could be merged with relatively little effort? I would agree with the OP that the (interval, real) method is the most common.

[Kolaru]

Since I am anyways only searching for roots for x bounded in [0, 2pi], would you say it is correct to simply remove the mod operation from the above function?

To me it looks like you could use mod at the very end of the computation, or not at all, since for what matters it goes into sin anyway. Feels like the mod2pi is here for numerical stability.

I merged #525 for now. It will not remove the error directly, because it doesn't add mod2pi. It can be avoided with mod(x, 2pi). Also it may cause a new problem : if you get an interval that overlap with one of the bound, the mod will get you the full range (e.g. mod(-0.0001 .. 0.0002, 2pi) == (0..2pi)),

What you want may rather be something like xn = mod2pi(mid(s)) ± radius(s), which doesn't guarantee inclusion in the range [0, 2pi], but will keep xn close enough, without risking a huge increase in the interval radius. Disclaimer : I haven't check if the formula correctly include the numerical error, I would need a bit more time to come up with something that guarantee it.

[Datseris]

Thank you yes. In the meantime I indeed went ahead and removed the mod from my equations but unfortunately in the end of the day it didn't help me, turns out the Jacobian matrix isn't invertible so it was all for naught xD But thanks a lot for progressing with the mod!