julia> mod(1.0, 0.1)
0.09999999999999995
Shouldn’t this be exactly zero, or at least pretty close to zero?
julia> mod(1.1, 0.1)
2.7755575615628914e-17
Like this, for instance?
The literal 0.1 isn’t 0.1, but the closest Float64 to it:
julia> using Printf
julia> @printf "%.60f" 0.1
0.100000000000000005551115123125782702118158340454101562500000
It is larger than 0.1, so
julia> fld(1.0,0.1)
9.0
julia> fma(-9.0, 0.1, 1.0) # 1.0 - 9.0 * 0.1 without rounding
0.09999999999999995
and yet, in Matlab:
>> mod(1.0, 0.1)
ans =
0
Matlab is lying to you.
4 Likes
However, something like this might surprise a Julia user:
julia> test(x,y) = x - y*round(x/y, RoundToZero) - rem(x,y)
test (generic function with 1 method)
julia> test(1.0, 0.1)
-0.09999999999999995
reading the documentation (even if it states without any intermediate rounding, which is not an obvious disclamer).
EDIT:
I am aware of:
julia> round(big(1.0)/big(0.1), RoundToZero)
9.0
and
julia> Float64(rem(Rational{BigInt}(1.0), Rational{BigInt}(0.1)))
0.09999999999999995
But this does not change that I think users might be confused, so maybe we should explain it better.
These are the tests to consider:
julia> test2(x,y) = x - fma(y, div(x, y), rem(x,y))
julia> test2(1.0,0.1)
0.0
julia> test3(x,y) = x - fma(y, fld(x, y), mod(x,y))
julia> test3(1.0,0.1)
0.0
I agree, but currently the documentation of rem uses round (that is why I have written that maybe we should improve it).
Neither of these tests does what I expect:
julia> x = 1.0
julia> y = 0.3
julia> x - fma(y, fld(x, y), mod(x,y))
0.0
To be clear, I’m looking for a function that returns zero when the value x is evenly divisible by y, and non-zero otherwise. Or, following the above point, when x is very close to being evenly divisible by y.
It appears that perhaps what I want is captured by:
x - trunc(x/y)*y
You will still have a problem if you take e.g. x = prevfloat(1.0) or y = nextfloat(0.1).
If you really want something that is very close to zero then maybe something like this will do what you want:
f(x,y) = min(abs(mod(x,y)), abs(mod(x,-y)))
which gives you a minimum deviation to being closely divisible from below and above.
It isn’t quite lying, it’s just being Matlab. From their docs,
Note that
modattempts to compensate for floating-point round-off effects to produce exact integer results when possible.
2 Likes
this is substantively faster, if it works for you: g(x,y) = mod(-abs(x), abs(y))
However, it is incorrect for e.g. g(1.0,1/3)
Also consider h(x,y) = (z=x/y; z-trunc(z) == 0.0 ? 0.0 : mod(x,y))
h(1.0, 0.1) yields 0.0 and h(1.0, 1/3) yields 0.0 too.
It seems to work, but do your own tests to verify.
Well, that depends how you define “evenly divisible” in the context of floating point. Perhaps if you explained the context in which you want it?
The problem with this definition is:
julia> g(1.0, 1/3)
0.33333333333333326
if I understand the original question correctly.
My specific problem is that I have implemented a model predictive control algorithm. This can be thought of as an algorithm that runs at a fast, hopefully-but-not-really-constant rate and calculates some output that gets sent to a hardware device. I want to output the output of the routine at some slower rate, but at a hopefully-nearly-constant time interval. As an example, my routine runs at approximately 100 iterations per second, and I want output at 10 iterations per second. Note that the output rate can vary at the user’s discretion, but should always be approximately an even divisor of the inner loop running rate.
Internally, elapsed time is available to me. The tempting thing to do, therefore, is to call time “t” and desired output rate “dt”, and call a routine to record the output whenever mod(t, dt) is approximately zero.
Also note that although I am using Julia, this is because I am the algorithm prototyper. The eventual target code will be implemented in C, by a hard-realtime software engineer who frowns on garbage collection, dynamic memory allocation, and really any advanced language features. So I am somewhat limited in the language features I can use. Iterators, for instance, are probably not viable here.
The other obvious way to do this, of course, would be to set up a counter and output when the counter reaches zero. This has the disadvantage that if my running rate ever varies (eg it slows down), my output rate will vary also. Which is not acceptable in my application.
1 Like
well – if you want to get technical 
here is one, about 2x the original with the advantage of giving 0.0 rather than close to zero for a number of inputs
function h(x,y)
a = x*inv(y)
a - trunc(a)
end
1 Like
@GlenHenshaw please watch the Julia Robotics talk at JuliaCom2018
They have the very same requirements you have - and a control loop which needs to run at 1000Hz
they had to avoid allocations in the inner tight loop.
Please listen to the questions at the end - where Valentin is asked about a macro which will avoid allocations.
And it is promised! 38:40 into the presentation, during questions
That seems dicey to me. If by “approximately zero” you mean abs(mod(t, dt)) < ε for some ε, then if you make ε too small you risk outputting very infrequently (because your code only gets called at a discrete set of times t, which may not often be close to multiples of dt), and if you make ε too large then you risk outputting multiple times in the same dt time interval.
Why not do something like:
t = time()
if t >= last_output_time + dt
...output...
last_output_time = t
end
That way it will call the output routine as close as possible to intervals of dt, limited by the frequency of the times t where you perform this check.
3 Likes
It might be possible to use the Rational type
julia> 1 % 1//10 == 0
true
1 Like
Very good point. Operations like these are actually one of the most important use cases of Rational. It is not in the Base by accident, or purely as a demo for parametric constructors 