Skip to content

Instantly share code, notes, and snippets.

@peerreynders

peerreynders/demo.exs

Last active Sep 8, 2019
Embed
What would you like to do?
Vampire Numbers Elixir demonstration script
# Cobbled together ("untested") demonstration code - no production value
# The intention is to flood all CPU cores with "work"
#
# Created for: https://elixirforum.com/t/finding-vampire-numbers/25098/10
#
defmodule Primes do
require Integer
@primes :primes
@init_last 3
@next_timeout 5_000
# ---
# GenServer related
def init(init_arg) do
precache = Keyword.get(init_arg, :precache, @init_last)
{:ok, [], {:continue, precache}}
end
def handle_continue(precache, []) do
initialize()
last_prime = :ets.last(@primes)
{:noreply, {%{}, last_prime, prime_next(precache, last_prime)}}
end
# keep generating values until limit is reached
#
def handle_cast(:next_prime, {requests, last_prime, precache}) do
new_prime = next_prime(last_prime)
new_requests = dispatch_replies(requests, new_prime)
new_precache = prime_next(precache, new_prime)
{:noreply, {new_requests, new_prime, new_precache}}
end
def handle_call({:next_prime, n}, from, {requests, last_prime, precache}) do
new_requests = add_request(requests, from, n)
new_precache =
case precache <= n do
false ->
# still to be generated
precache
_ when precache <= last_prime ->
# restart generation
prime_next(n + 1, last_prime)
_ ->
# bump up ongoing generation limit
n + 1
end
{:noreply, {new_requests, last_prime, new_precache}}
end
def start_link(opts) do
{precache, opts} = Keyword.pop(opts, :precache, @init_last)
init_arg = [precache: precache]
opts = [{:name, __MODULE__} | opts]
GenServer.start_link(__MODULE__, init_arg, opts)
end
# ---
# GenServer State
defp prime_next(precache, last_prime) when precache > last_prime do
GenServer.cast(self(), :next_prime)
precache
end
defp prime_next(_, last_prime) do
IO.puts("Primes stopped: #{last_prime}")
last_prime
end
defp add_request(requests, from, last_prime),
do: Map.update(requests, last_prime, [from], &[from | &1])
defp dispatch_replies(requests, last_prime) do
requests
|> Map.to_list()
|> List.foldl(%{}, make_dispatch(last_prime))
end
defp make_dispatch(last_prime) do
fn
{n, pids}, pending when n < last_prime ->
next_prime = :ets.next(@primes, n)
dispatch_next(next_prime, pids)
pending
{n, pids}, pending ->
# not ready
Map.put(pending, n, pids)
end
end
defp dispatch_next(_, []) do
:ok
end
defp dispatch_next(value, [pid | others]) do
GenServer.reply(pid, value)
dispatch_next(value, others)
end
# ---
# Cache population
defp initialize() do
:ets.new(@primes, [:ordered_set, :named_table])
:ets.insert(@primes, {2})
:ets.insert(@primes, {@init_last})
end
defp next_prime(n) when Integer.is_odd(n),
do: add_next_prime(n + 2)
defp add_next_prime(n) do
if uncached_prime?(n) do
:ets.insert(@primes, {n})
n
else
add_next_prime(n + 2)
end
end
defp uncached_prime?(n),
do: uncached_prime?(n, @init_last)
# number divisible by a prime cannot be a prime
defp uncached_prime?(n, p) when rem(n, p) == 0 do
false
end
defp uncached_prime?(n, last_prime) do
case :ets.next(@primes, last_prime) do
prime when is_integer(prime) ->
uncached_prime?(n, prime)
_ ->
# ran out of primes - n as prime
true
end
end
# ---
# Client functions
def first(),
do: :ets.first(@primes)
def next(last_prime) do
case :ets.next(@primes, last_prime) do
prime when is_integer(prime) ->
prime
_ ->
prime_request(last_prime)
end
end
defp prime_request(last_prime),
do: GenServer.call(__MODULE__, {:next_prime, last_prime}, @next_timeout)
end
defmodule Vnum do
require Integer
def isqrt(x) when x < 0,
do: raise(ArithmeticError)
def isqrt(x),
do: isqrt(x, 1, div(1 + x, 2))
defp isqrt(x, m, n) when abs(m - n) <= 1 and n * n <= x,
do: n
defp isqrt(_x, m, n) when abs(m - n) <= 1,
do: n - 1
defp isqrt(x, _, n),
do: isqrt(x, n, div(n + div(x, n), 2))
#
# prime_factors(1530, Primes, 2, [])
#
# Output:
# [{2, 1}, {3, 2}, {5, 1}, {17, 1}]
#
defp prime_factors(num, primes, prime, factors) do
if prime <= isqrt(num) do
{exponent, remainder} = prime_exponent(num, prime, 0)
new_factors = if(exponent > 0, do: [{prime, exponent} | factors], else: factors)
next_prime = primes.next(prime)
prime_factors(remainder, primes, next_prime, new_factors)
else
new_factors = if(num >= 2, do: [{num, 1} | factors], else: factors)
:lists.reverse(new_factors)
end
end
defp prime_exponent(num, prime, count) when rem(num, prime) == 0 do
num
|> div(prime)
|> prime_exponent(prime, count + 1)
end
defp prime_exponent(num, _, count) do
{count, num}
end
#
# Input:
# [{2, 1}, {3, 2}, {5, 1}, {17, 1}]
#
# Output:
# [17, 85, 5, 153, 765, 45, 9, 51, 255, 15, 3, 34, 170, 10, 306, 1530, 90, 18, 102, 510, 30, 6, 2, 1]
#
defp divisors(factors),
do: divisors(factors, 1, [1])
defp divisors([], _, acc),
do: acc
defp divisors([{prime, exponent} | rest], factor, acc) do
new_acc = acc_divisors(rest, prime, exponent, factor, acc)
divisors(rest, factor, new_acc)
end
defp acc_divisors(factors, prime, count, factor, acc) when count > 0 do
new_factor = factor * prime
new_acc = divisors(factors, new_factor, [new_factor | acc])
acc_divisors(factors, prime, count - 1, new_factor, new_acc)
end
defp acc_divisors(_, _, _, _, acc) do
acc
end
#
# Input:
# [17, 85, 5, 153, 765, 45, 9, 51, 255, 15, 3, 34, 170, 10, 306, 1530, 90, 18, 102, 510, 30, 6, 2, 1]
#
# Output:
# [{1, 1530}, {2, 765}, {3, 510}, {5, 306}, {6, 255}, {9, 170}, {10, 153}, {15, 102}, {17, 90}, {18, 85}, {30, 51}, {34, 45}]
#
defp pair_divisors(divisors) do
count =
divisors
|> length()
|> div(2)
{lower, upper} =
divisors
|> Enum.sort()
|> Enum.split(count)
List.zip([lower, :lists.reverse(upper)])
end
#
# ordered_digits(1530) -> {"0135", 4}
#
defp ordered_digits(n) do
digits = Integer.digits(n)
{ordered_string(digits), length(digits)}
end
#
# ordered_string([1, 5, 3, 0]) -> "0135"
#
defp ordered_string(digits) do
digits
|> Enum.sort()
|> Enum.join()
end
# ref: https://oeis.org/A014575
#
# fangs?("0135", 2, {30, 51}) -> true
# fangs?("0135", 2, {34, 45}) -> false
#
defp fangs?(ordered, size, {x, y}) do
with x_digits <- Integer.digits(x),
{x_digits, true} <- {x_digits, length(x_digits) === size},
{x_digits, y_digits} <- {x_digits, Integer.digits(y)},
{x_digits, y_digits, true} <- {x_digits, y_digits, length(y_digits) === size},
{x_digits, y_digits, true} <- {x_digits, y_digits, rem(x, 10) !== 0 or rem(y, 10) !== 0} do
ordered == ordered_string(x_digits ++ y_digits)
else
_ ->
false
end
end
def filter_fangs(pairs, ordered, fang_size),
do: Enum.filter(pairs, &fangs?(ordered, fang_size, &1))
def extract_fangs(primes, n) do
{ordered, digit_count} = ordered_digits(n)
with true <- Integer.is_even(digit_count) and digit_count > 2 do
fang_size = div(digit_count, 2)
n
|> prime_factors(primes, 2, [])
|> divisors()
|> pair_divisors()
|> filter_fangs(ordered, fang_size)
else
_ ->
[]
end
end
end
defmodule Demo do
def run(x, y) do
{m, n} = adjust(x, y)
Primes.start_link(precache: div(n, 2))
opts = [timeout: 5_000]
m..n
|> Task.async_stream(__MODULE__, :extract_fangs, [Primes], opts)
|> Enum.reduce([], &keep_vampires/2)
end
defp adjust(x, y) when is_integer(x) and is_integer(y) do
{m, n} = if(x <= y, do: {x, y}, else: {y, x})
{max(m, 1000), max(n, 1000)}
end
def extract_fangs(n, primes),
do: {n, Vnum.extract_fangs(primes, n)}
defp keep_vampires({:ok, {_n, [_ | _]} = vampire_result}, others),
do: [vampire_result | others]
defp keep_vampires(_non_vampire_result, others),
do: others
end
IO.puts("Schedulers online: #{System.schedulers_online()}")
# http://oeis.org/A048939
args = [1_000_100_000, 1_001_300_000]
# args = [1_000_000_000, 2_000_000_000]
{time, result} = :timer.tc(Demo, :run, args)
IO.puts("Seconds: #{div(time, 1_000_000)}")
IO.inspect(result)
$ elixir demo.exs
Schedulers online: 8
Seconds: 265
[
  {1001295697, [{19001, 52697}]},
  {1001288794, [{10874, 92081}]},
  {1001249680, [{10421, 96080}]},
  {1001244420, [{24420, 41001}]},
  {1001240352, [{25011, 40032}]},
  {1001193790, [{10310, 97109}]},
  {1001147422, [{24002, 41711}]},
  {1001139282, [{11022, 90831}]},
  {1001127442, [{14021, 71402}]},
  {1001020252, [{20012, 50021}]},
  {1000812510, [{12510, 80001}]},
  {1000520064, [{20004, 50016}]},
  {1000520010, [{20010, 50001}]},
  {1000425010, [{25010, 40001}]},
  {1000407528, [{12504, 80007}]},
  {1000250010, [{20001, 50010}]},
  {1000198206, [{10206, 98001}]},
  {1000191991, [{10991, 91001}]},
  {1000174288, [{14288, 70001}]}
]
$ 
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.