Implement Euler totient function (#43)

Author: TotalVerb

src/Primes.jl

@@ -14,7 +14,7 @@ end
 using Base: BitSigned
 using Base.Checked.checked_neg
 
-export ismersenneprime, isrieselprime, nextprime, prevprime, prime, prodfactors, radical
+export ismersenneprime, isrieselprime, nextprime, prevprime, prime, prodfactors, radical, totient
 
 include("factorization.jl")
 
@@ -500,6 +500,37 @@ function ll_primecheck(X::Integer, s::Integer = 4)
     return S == 0
 end
 
+"""
+    totient(f::Factorization{T}) -> T
+
+Compute the Euler totient function of the number whose prime factorization is
+given by `f`. This method may be preferable to [`totient(::Integer)`](@ref)
+when the factorization can be reused for other purposes.
+"""
+function totient{T <: Integer}(f::Factorization{T})
+    if !isempty(f) && first(first(f)) == 0
+        throw(ArgumentError("ϕ(0) is not defined"))
+    end
+    result = one(T)
+    for (p, k) in f
+        result *= p^(k-1) * (p - 1)
+    end
+    result
+end
+
+"""
+    totient(n::Integer) -> Integer
+
+Compute the Euler totient function ``ϕ(n)``, which counts the number of
+positive integers less than or equal to ``n`` that are relatively prime to
+``n`` (that is, the number of positive integers `m ≤ n` with `gcd(m, n) == 1`).
+The totient function of `n` when `n` is negative is defined to be
+`totient(abs(n))`.
+"""
+function totient(n::Integer)
+    totient(factor(abs(n)))
+end
+
 # add_! : "may" mutate the Integer argument (only for BigInt currently)
 
 # modify a BigInt in-place

test/runtests.jl

@@ -276,6 +276,37 @@ for T = (Int, UInt, BigInt)
     @test issorted(f.pe)
 end
 
+# dumb implementation of Euler totient for correctness tests
+ϕ(n) = count(m -> gcd(n, m) == 1, 1:n)
+
+@testset "totient(::$T)" for T in [Int16, Int32, Int64, BigInt]
+    for n in 1:1000
+        @test ϕ(T(n)) == totient(T(n)) == totient(-T(n))
+    end
+    @test_throws ArgumentError @inferred(totient(T(0)))
+end
+
+# check some big values
+@testset "totient() correctness" begin
+    @test totient(2^4 * 3^4 * 5^4) == 216000
+    @test totient(big"2"^1000) == big"2"^999
+
+    const some_coprime_numbers = BigInt[
+        450000000, 1099427429702334733, 200252151854851, 1416976291499, 7504637909,
+        1368701327204614490999, 662333585807659, 340557446329, 1009091
+    ]
+
+    for i in some_coprime_numbers
+        for j in some_coprime_numbers
+            if i ≠ j
+                @test totient(i*j) == totient(i) * totient(j)
+            end
+        end
+        # can use directly with Factorization
+        @test totient(i) == totient(factor(i))
+    end
+end
+
 # check copy property for big primes relied upon in nextprime/prevprime
 for n = rand(big(-10):big(10), 10)
     @test n+0 !== n