I want to compute the integer m in m = log(b, y) for b integer, y=b^m.
It would return an error if y is not a power of b.
Right now, using the log is not ideal since
julia> log(7,7^5)
5.000000000000001
which is not a Int.
I see ways to code a function f(b,y) doing that, I just wondered if there were a dedicated function already doing that.
If you know that the result for your arguments ought to be an integer, just define an integer log function like this
julia> logi = Int∘round∘log
Int64 ∘ round ∘ log
julia> logi(7,7^5)
5
or for better type safety:
julia> logi(b::Int, y::Int) = Int(round(log(b,y)))
logi (generic function with 1 method)
julia> logi(3.5,100.0)
ERROR: MethodError: no method matching logi(::Float64, ::Float64)
Stacktrace:
[1] top-level scope
@ REPL[3]:1
julia> logi(7,7^5)
5
You could even do type piracy if it’s for your own use only:
julia> Base.log(b::Int, y::Int) = Int(round(log(Float64(b),Float64(y))))
julia> log(7,7^5)
5
imo we should have this in base because it’s easy to implement badly and hard to do well.
Does the Float as an intermediate is the best strategy for performance? ps: Maybe a logi function in base doing that should throw an error if the result is not exactly representable by an integer.
ps: @Jeff_Emanuel, why not
logi(x::Int,y::Int) = round(Int, log(x,y))
?
Just taste, or is there any other reason behind?
Now the problem is that
julia> logi(7,8^5)
5
where i would expect an error.
I share this feeling
No particular reason. It’s just more natural for me as Int(round(... rather than round(Int,....
Certainly. I thought from your description that in your use case you knew that the second argument was an integral power of your base, thus the overly simplistic suggestion.
I want to both compute m and check that it is an integer.
If log(b,y) was exact it would return exactly m when y = b^m so I would just do m = Int(log(b,y)).
However, it seems that in some case this do not work e.g. m = Int(log(b,7,7^5)).
That is trickier than it seems.
For instance:
julia> convert(Int, log(7,7^5))
ERROR: InexactError: Int64(5.000000000000001)
Thus, the floating point log breaks the representability of the result as an Integer.
The function should probably be something on the lines of:
julia> function logi(x::Int,y::Int)
x == 0 && error()
i = 0
r = x
while abs(r) < abs(y)
i += 1
r *= x
end
r == y ? i+1 : error()
end
logi (generic function with 2 methods)
julia> logi(-7,(-7)^5)
5
julia> logi(7,(7)^5)
5
but written more carefully, because that one fails in a series of cases, for instance:
julia> logi(-1,1)
ERROR:
I would definitely like a floorlog/ceillog function. Note that elements of these hypothetical functions already exist in the implementations of prevpow/nextpow.
Or maybe these are implemented as RoundingMode arguments to Base.log? But that might be piling too much disparate functionality under one function.
julia> function logi(b::Int, y::Int)
n = round(Int, log(b, y))
b^n == y || throw(ArgumentError("$y is not a power of $b"))
return n
end
logi (generic function with 1 method)
julia> using BenchmarkTools
julia> b = 5; y = 5^12;
julia> @btime logi($b, $y)
32.060 ns (0 allocations: 0 bytes)
12
For comparison:
julia> @btime log($b, $y)
35.000 ns (0 allocations: 0 bytes)
12.000000000000002
Edit: Perhaps someone could explain why logi is slightly but consistently faster than log with the same inputs.
Isn’t this a hard problem to solve exactly in general, i.e. the discrete logarithm problem if you want to support BigInt? I guess it only becomes hard if there is a mod operation?
Yeah. It’s only hard with a modulo. Without, there is a simple algorithm that’s O(log(n)) steps that just repeatedly squares until you pass the number and than does a binary search to find the true answer.