first commit

util/julia/files/_erastothenes_extended_sieve.jl

@@ -0,0 +1,59 @@
+include("_erastothenes_sieve.jl")
+include("_math.jl")
+
+mutable struct extended_sieve
+    sieve::sieve
+    factors::Vector{Vector{Int}}
+end
+
+function get_extended_sieve()
+    return extended_sieve(sieve(falses(1), 1, 1, 1000), [[]])
+end
+
+function _elongate!(es::extended_sieve)
+    s = es.sieve
+    length = s._length
+    target = s._target
+    diff = target - length
+    run!(s, target)
+
+    factors = append!(es.factors, fill([], nneg(diff)))
+    done = [trues(length); falses(diff)]
+    if diff > 0
+        for i ∈ get_primes(s)
+            if i > length
+                factors[i] = [i]
+            end
+            for j ∈ cld(length+1, i):fld(target, i)
+                prod = i*j
+                if ~done[prod]
+                    if done[j]
+                        index = searchsortedfirst(factors[j], i)
+                        factors[prod] = insert!(copy(factors[j]), index, i)
+                        done[prod] = true
+                    else
+                        push!(factors[prod], i)
+                    end
+                end
+            end
+        end
+    end
+    return Nothing
+end
+
+function run!(es::extended_sieve, n::Int)
+    es.sieve._target = n;
+    _elongate!(es)
+    return Nothing
+end
+
+function get(es::extended_sieve, n::Int, prerun=true::Bool)
+    if n > es.sieve._length
+        run(s, prerun ? n + es.sieve._prerun_size : n)
+    end
+    return (es.sieve[n], es.factors[n])
+end
+
+function get_primes(es::extended_sieve)
+    return get_primes(es.sieve)
+end

util/julia/files/_erastothenes_sieve.jl

@@ -0,0 +1,88 @@
+mutable struct sieve
+    sieve::BitVector
+    _length::Int
+    _target::Int
+    _prerun_size::Int
+end
+
+function get_sieve()
+    return sieve([false], 1, 1, 1000)
+end
+
+function _elongate!(s::sieve)
+    length = s._length
+    target = s._target
+    diff = target - length
+    if diff > 0
+        s.sieve = [s.sieve; trues(diff)]
+        for i = 1:floor(Int, sqrt(target))
+            if s.sieve[i]
+                _floor = length÷i + 1
+                _ceil = target÷i
+                for j = _floor:_ceil
+                    actual = i*j
+                    if actual > i
+                        s.sieve[actual] = false
+                    end
+                end
+            end
+        end
+        s._length = target
+    end
+    return nothing
+end
+
+# runs the sieve until index n
+function run!(s::sieve, n::Int)
+    s._target = n
+    _elongate!(s)
+    return nothing
+end
+
+# fetches data about a specific integer
+function get(s::sieve, n::Int, prerun=true::Bool)
+    if n > s._length
+        run(s, prerun ? n + s._prerun_size : n)
+    end
+    return s.sieve[n]
+end
+
+# returns all primes in the sieve
+function get_primes(s::sieve, n_max=(-1)::Int)
+    if n_max == -1
+        trimmed_sieve = s.sieve
+    else
+        trimmed_sieve = s.sieve[1:n_max]
+    end
+    n = length(trimmed_sieve)
+    c = count(trimmed_sieve)
+    primes = Vector(undef, c)
+    j = 1
+    for i = 1:n
+        if s.sieve[i]
+            primes[j] = i
+            j += 1
+        end
+    end
+    return primes
+end
+
+function get_composites(s::sieve)
+    n = length(s.sieve)
+    c = count(s.sieve)
+    composites = Vector(undef, n-c)
+    j = 1
+    for i = 1:n
+        if !s.sieve[i]
+            composites[j] = i
+            j += 1
+        end
+    end
+    return composites
+end
+
+function is_prime(s::sieve, i::Int)
+    if i < 0 return false end
+    run(s, i)
+    return s.sieve[i]
+end

util/julia/files/_math.jl

@@ -0,0 +1,3 @@
+
+nneg(v) = max(v, 0)
+npos(v) = min(v, 0)

util/julia/files/_sorting.jl

@@ -0,0 +1,31 @@
+
+function mergeSorted(A::Vector{T}, B::Vector{T}) where T
+    la = length(A)
+    lb = length(B)
+    C = Vector{T}(undef, la + lb)
+
+    a = b = 1
+    while (a <= la && b <= lb)
+        if A[a] < B[b]
+            C[a+b-1] = A[a]
+            a += 1
+        else
+            C[a+b-1] = B[b]
+            b += 1
+        end
+    end
+
+    if a > la
+        while b <= lb
+            C[a+b-1] = B[b]
+            b += 1
+        end
+    else
+        while a <= la
+            C[a+b-1] = A[a]
+            a += 1
+        end
+    end
+
+    return C
+end