From 175628c4d321e44af4adaffea677ca7da8266707 Mon Sep 17 00:00:00 2001 From: Rafael Fourquet Date: Tue, 3 Dec 2019 10:35:50 +0100 Subject: [PATCH 1/2] swith argument order in factor! and pollardfactors! This put the mutated argument in first position, as is recommended. --- src/Primes.jl | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/src/Primes.jl b/src/Primes.jl index 5b91d79..949abc4 100644 --- a/src/Primes.jl +++ b/src/Primes.jl @@ -235,7 +235,7 @@ isprime(n::Int128) = n < 2 ? false : # https://en.wikipedia.org/wiki/Pollard%27s_rho_algorithm # http://maths-people.anu.edu.au/~brent/pub/pub051.html # -function factor!(n::T, h::AbstractDict{K,Int}) where {T<:Integer,K<:Integer} +function factor!(h::AbstractDict{K,Int}, n::T) where {T<:Integer,K<:Integer} # check for special cases if n < 0 h[-1] = 1 @@ -243,7 +243,7 @@ function factor!(n::T, h::AbstractDict{K,Int}) where {T<:Integer,K<:Integer} h[2] = 8 * sizeof(T) - 1 return h else - return factor!(checked_neg(n), h) + return factor!(h, checked_neg(n)) end elseif n == 1 return h @@ -265,7 +265,7 @@ function factor!(n::T, h::AbstractDict{K,Int}) where {T<:Integer,K<:Integer} isprime(n) && (h[n] = 1; return h) end end - T <: BigInt || widemul(n - 1, n - 1) ≤ typemax(n) ? pollardfactors!(n, h) : pollardfactors!(widen(n), h) + T <: BigInt || widemul(n - 1, n - 1) ≤ typemax(n) ? pollardfactors!(h, n) : pollardfactors!(h, widen(n)) end @@ -298,7 +298,7 @@ julia> collect(factor(0)) 0=>1 ``` """ -factor(n::T) where {T<:Integer} = factor!(n, Factorization{T}()) +factor(n::T) where {T<:Integer} = factor!(Factorization{T}(), n) """ @@ -337,7 +337,7 @@ julia> factor(Set, 100) Set([2,5]) ``` """ -factor(::Type{D}, n::T) where {T<:Integer, D<:AbstractDict} = factor!(n, D(Dict{T,Int}())) +factor(::Type{D}, n::T) where {T<:Integer, D<:AbstractDict} = factor!(D(Dict{T,Int}()), n) factor(::Type{A}, n::T) where {T<:Integer, A<:AbstractArray} = A(factor(Vector{T}, n)) factor(::Type{Vector{T}}, n::T) where {T<:Integer} = mapreduce(collect, vcat, [repeated(k, v) for (k, v) in factor(n)], init=Vector{T}()) @@ -383,7 +383,7 @@ julia> radical(2*2*3) """ radical(n) = prod(factor(Set, n)) -function pollardfactors!(n::T, h::AbstractDict{K,Int}) where {T<:Integer,K<:Integer} +function pollardfactors!(h::AbstractDict{K,Int}, n::T) where {T<:Integer,K<:Integer} while true c::T = rand(1:(n - 1)) G::T = 1 @@ -423,9 +423,9 @@ function pollardfactors!(n::T, h::AbstractDict{K,Int}) where {T<:Integer,K<:Inte G = gcd(x > ys ? x - ys : ys - x, n) end if G != n - isprime(G) ? h[G] = get(h, G, 0) + 1 : pollardfactors!(G, h) + isprime(G) ? h[G] = get(h, G, 0) + 1 : pollardfactors!(h, G) G2 = div(n,G) - isprime(G2) ? h[G2] = get(h, G2, 0) + 1 : pollardfactors!(G2, h) + isprime(G2) ? h[G2] = get(h, G2, 0) + 1 : pollardfactors!(h, G2) return h end end From 21e380a60ead1aa249231c75c784b055db075ee8 Mon Sep 17 00:00:00 2001 From: Rafael Fourquet Date: Tue, 3 Dec 2019 11:20:49 +0100 Subject: [PATCH 2/2] extend factor! to other containers and export it --- docs/src/api.md | 1 + src/Primes.jl | 77 ++++++++++++++++++++++++++++++++++++++---------- test/runtests.jl | 6 ++++ 3 files changed, 68 insertions(+), 16 deletions(-) diff --git a/docs/src/api.md b/docs/src/api.md index 5e08842..54174b5 100644 --- a/docs/src/api.md +++ b/docs/src/api.md @@ -8,6 +8,7 @@ DocTestSetup = :(using Primes) ```@docs Primes.factor +Primes.factor! Primes.prodfactors ``` diff --git a/src/Primes.jl b/src/Primes.jl index 949abc4..0f62a9b 100644 --- a/src/Primes.jl +++ b/src/Primes.jl @@ -8,7 +8,7 @@ import Base: iterate, eltype, IteratorSize, IteratorEltype using Base: BitSigned using Base.Checked: checked_neg -export isprime, primes, primesmask, factor, ismersenneprime, isrieselprime, +export isprime, primes, primesmask, factor, factor!, ismersenneprime, isrieselprime, nextprime, nextprimes, prevprime, prevprimes, prime, prodfactors, radical, totient include("factorization.jl") @@ -229,43 +229,88 @@ isprime(n::Int128) = n < 2 ? false : n ≤ typemax(Int64) ? isprime(Int64(n)) : isprime(BigInt(n)) +push_factor!(h::AbstractVector, f, mult) = append!(h, Iterators.repeated(f, mult)) +push_factor!(h::AbstractVector, f) = push!(h, f) +push_factor!(h::AbstractDict, f, mult=1) = h[f] = get(h, f, 0) + mult +push_factor!(h::AbstractSet, f, mult=1) = push!(h, f) + +maybe_sort!(h) = h +maybe_sort!(h::AbstractVector) = sort!(h) + +FactorCont{K} = Union{AbstractVector{K}, AbstractDict{K,Int}, AbstractSet{K}} where K + # Trial division of small (< 2^16) precomputed primes + # Pollard rho's algorithm with Richard P. Brent optimizations # https://en.wikipedia.org/wiki/Trial_division # https://en.wikipedia.org/wiki/Pollard%27s_rho_algorithm # http://maths-people.anu.edu.au/~brent/pub/pub051.html # -function factor!(h::AbstractDict{K,Int}, n::T) where {T<:Integer,K<:Integer} +""" + factor!(container, n::Integer) -> container + +Return the prime factorization of `n` stored in `container`. +When `container::AbstractDict`, the keys represent the factors and the values +represent the multiplicities. +When `container::AbstractSet`, the elements represent the distinct factors +of `n` (the multiplicities are not stored). +When `container::AbstractVector`, the factors are stored with their multiplicities, +in sorted order. + +# Examples +```julia +julia> factor!(Dict{Int,Int}(), 63) +Dict{Int64,Int64} with 2 entries: + 7 => 1 + 3 => 2 + +julia> factor!(Set{Int}(), 63) +Set{Int64} with 2 elements: + 7 + 3 + +julia> factor!(Int[], 63) +3-element Array{Int64,1}: + 3 + 3 + 7 + +julia> prod(ans) +63 +``` +""" +function factor!(h::FactorCont{K}, n::T) where {T<:Integer,K<:Integer} # check for special cases if n < 0 - h[-1] = 1 + push_factor!(h, -1) if isa(n, BitSigned) && n == typemin(T) - h[2] = 8 * sizeof(T) - 1 - return h + push_factor!(h, 2, 8 * sizeof(T) - 1) + return maybe_sort!(h) else return factor!(h, checked_neg(n)) end elseif n == 1 - return h + return maybe_sort!(h) elseif n == 0 || isprime(n) - h[n] = 1 - return h + push_factor!(h, n) + return maybe_sort!(h) end local p::T for p in PRIMES if n % p == 0 - h[p] = get(h, p, 0) + 1 + push_factor!(h, p) n = div(n, p) while n % p == 0 - h[p] = get(h, p, 0) + 1 + push_factor!(h, p) n = div(n, p) end - n == 1 && return h - isprime(n) && (h[n] = 1; return h) + n == 1 && return maybe_sort!(h) + isprime(n) && (push_factor!(h, n); return h) end end - T <: BigInt || widemul(n - 1, n - 1) ≤ typemax(n) ? pollardfactors!(h, n) : pollardfactors!(h, widen(n)) + maybe_sort!(T <: BigInt || widemul(n - 1, n - 1) ≤ typemax(n) ? + pollardfactors!(h, n) : + pollardfactors!(h, widen(n))) end @@ -383,7 +428,7 @@ julia> radical(2*2*3) """ radical(n) = prod(factor(Set, n)) -function pollardfactors!(h::AbstractDict{K,Int}, n::T) where {T<:Integer,K<:Integer} +function pollardfactors!(h::FactorCont{K}, n::T) where {T<:Integer,K<:Integer} while true c::T = rand(1:(n - 1)) G::T = 1 @@ -423,9 +468,9 @@ function pollardfactors!(h::AbstractDict{K,Int}, n::T) where {T<:Integer,K<:Inte G = gcd(x > ys ? x - ys : ys - x, n) end if G != n - isprime(G) ? h[G] = get(h, G, 0) + 1 : pollardfactors!(h, G) + isprime(G) ? push_factor!(h, G) : pollardfactors!(h, G) G2 = div(n,G) - isprime(G2) ? h[G2] = get(h, G2, 0) + 1 : pollardfactors!(h, G2) + isprime(G2) ? push_factor!(h, G2) : pollardfactors!(h, G2) return h end end diff --git a/test/runtests.jl b/test/runtests.jl index 80728bc..f3a556f 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -239,6 +239,12 @@ end # factor returns a sorted dict @test all([issorted(collect(factor(rand(Int)))) for x in 1:100]) +@testset "factor!" begin + @test factor!(Int[], -50) == [-1, 2, 5, 5] + @test factor!(Set{Int}(), -50) == Set([-1, 2, 5]) + @test factor!(Dict{Int,Int}(), -50) == Dict(-1 => 1, 2 => 1, 5 => 2) +end + # Lucas-Lehmer @test !ismersenneprime(2047) @test ismersenneprime(8191)