Eisenstein primes

From Rosetta Code
Eisenstein primes is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

An Eisenstein integer is a non-unit Gaussian integer a + bω where ω(1+ω) = -1, and a and b are integers.

ω is generally chosen as a cube root of unity:

As a Gaussian integer, any Eisenstein integer is either a unit (an integer with a multiplicative inverse [±1, ±ω, ±(ω^-1)]), a prime (a number p such that if p divides xy, then p necessarily divides either x or y), or composite (a product of primes).

An Eisenstein integer a + bω is a prime if either it is a product of a unit and an integer prime p such that p % 3 == 2 or norm(a + bω) is an integer prime.

Eisenstein numbers can be generated by choosing any a and b such that a and b are integers. To allow generation in a relatively fixed order, such numbers can be ordered by their 2-norm or norm:

    norm(eisenstein integer(a, b)) = |a + bω| = 
Task
  • Find, and show here as their complex number values, the first 100 (by norm order) Eisenstein primes nearest 0.
Stretch
  • Plot in the complex plane at least the first 2000 such numbers (again, as found with norm closest to 0).
See also


Julia

""" rosettacode.org/wiki/Eisenstein_primes """

import Base: Complex, real, imag
import LinearAlgebra: norm
import Primes: isprime
import Plots: scatter

struct Eisenstein{T<:Integer} <: Number
    a::T
    b::T
    Eisenstein(a::T, b::T) where {T} = new{T}(a, b)
    Eisenstein(a::T) where {T<:Real} = new{T}(a, zero(T))
    Eisenstein(a::Real, b::Real) = new{eltype(promote(a, b))}(promote(a, b)...)
end

const ω = Eisenstein(false, true)

real(n::Eisenstein) = n.a - n.b / 2

imag(n::Eisenstein) = n.b * sqrt(big(3.0)) / 2

norm(n::Eisenstein) = n.a * n.a + n.b * n.b - n.a * n.b

Complex(n::Eisenstein{T}) where {T} = Complex{typeof(real(n))}(real(n), imag(n))

"""
    is_eisenstein_prime(n)

An Eisenstein integer is a non-unit Gaussian integer a + bω where ω(1+ω) = -1,
and a and b are integers. As a Gaussian integer, any Eisenstein integer is
either a unit (an integer with a multiplicative inverse [±1, ±ω, ±(ω^-1)]),
prime (a number p such that if p divides xy, then p necessarily divides
either x or y), or composite (a product of primes).

    An Eisenstein integer a + bω is a prime if either it is a product of a unit
and an integer prime p such that p % 3 == 2 or norm(a + bω) is an integer prime.
"""
function is_eisenstein_prime(n::Eisenstein)
    if n.a == 0 || n.b == 0 || n.a === n.b
        c = max(abs(n.a), abs(n.b))
        return isprime(c) && c % 3 == 2
    else
        return isprime(norm(n))
    end
end

function test_eisenstein_primes(graphlimitsquared = 10_000, printlimit = 100)
    lim = isqrt(graphlimitsquared)
    arr = [Eisenstein(a, b) for a = -lim:lim, b = -lim:lim]
    eprimes = sort!([arr[p[1]] for p in enumerate(arr) if is_eisenstein_prime(p[2])],
       lt = (x, y) -> norm(x) < norm(y))
    for (i, c) in enumerate(eprimes)
        if i <= printlimit
            print(lpad(Complex{Float32}(Complex(c)), 24), i % 4 == 0 ? "\n" : "")
        end
    end
    display(
        scatter(
            map(real, eprimes),
            map(imag, eprimes),
            markersize = 1,
            title = "Eisenstein primes with norm < $lim",
        ),
    )
end

test_eisenstein_primes()
Output:
   0.0f0 - 1.7320508f0im  -1.5f0 - 0.8660254f0im   1.5f0 - 0.8660254f0im  -1.5f0 + 0.8660254f0im
   1.5f0 + 0.8660254f0im   0.0f0 + 1.7320508f0im  -1.0f0 - 1.7320508f0im   1.0f0 - 1.7320508f0im
        -2.0f0 + 0.0f0im         2.0f0 + 0.0f0im  -1.0f0 + 1.7320508f0im   1.0f0 + 1.7320508f0im
   -0.5f0 - 2.598076f0im    0.5f0 - 2.598076f0im  -2.0f0 - 1.7320508f0im   2.0f0 - 1.7320508f0im
  -2.5f0 - 0.8660254f0im   2.5f0 - 0.8660254f0im  -2.5f0 + 0.8660254f0im   2.5f0 + 0.8660254f0im
  -2.0f0 + 1.7320508f0im   2.0f0 + 1.7320508f0im   -0.5f0 + 2.598076f0im    0.5f0 + 2.598076f0im
  -1.0f0 - 3.4641016f0im   1.0f0 - 3.4641016f0im   -2.5f0 - 2.598076f0im    2.5f0 - 2.598076f0im
  -3.5f0 - 0.8660254f0im   3.5f0 - 0.8660254f0im  -3.5f0 + 0.8660254f0im   3.5f0 + 0.8660254f0im
   -2.5f0 + 2.598076f0im    2.5f0 + 2.598076f0im  -1.0f0 + 3.4641016f0im   1.0f0 + 3.4641016f0im
  -0.5f0 - 4.3301272f0im   0.5f0 - 4.3301272f0im   -3.5f0 - 2.598076f0im    3.5f0 - 2.598076f0im
  -4.0f0 - 1.7320508f0im   4.0f0 - 1.7320508f0im  -4.0f0 + 1.7320508f0im   4.0f0 + 1.7320508f0im
   -3.5f0 + 2.598076f0im    3.5f0 + 2.598076f0im  -0.5f0 + 4.3301272f0im   0.5f0 + 4.3301272f0im
  -2.5f0 - 4.3301272f0im   2.5f0 - 4.3301272f0im        -5.0f0 + 0.0f0im         5.0f0 + 0.0f0im
  -2.5f0 + 4.3301272f0im   2.5f0 + 4.3301272f0im   -2.0f0 - 5.196152f0im    2.0f0 - 5.196152f0im
  -3.5f0 - 4.3301272f0im   3.5f0 - 4.3301272f0im  -5.5f0 - 0.8660254f0im   5.5f0 - 0.8660254f0im
  -5.5f0 + 0.8660254f0im   5.5f0 + 0.8660254f0im  -3.5f0 + 4.3301272f0im   3.5f0 + 4.3301272f0im
   -2.0f0 + 5.196152f0im    2.0f0 + 5.196152f0im  -0.5f0 - 6.0621777f0im   0.5f0 - 6.0621777f0im
  -5.0f0 - 3.4641016f0im   5.0f0 - 3.4641016f0im   -5.5f0 - 2.598076f0im    5.5f0 - 2.598076f0im
   -5.5f0 + 2.598076f0im    5.5f0 + 2.598076f0im  -5.0f0 + 3.4641016f0im   5.0f0 + 3.4641016f0im
  -0.5f0 + 6.0621777f0im   0.5f0 + 6.0621777f0im  -2.5f0 - 6.0621777f0im   2.5f0 - 6.0621777f0im
   -4.0f0 - 5.196152f0im    4.0f0 - 5.196152f0im  -6.5f0 - 0.8660254f0im   6.5f0 - 0.8660254f0im
  -6.5f0 + 0.8660254f0im   6.5f0 + 0.8660254f0im   -4.0f0 + 5.196152f0im    4.0f0 + 5.196152f0im
  -2.5f0 + 6.0621777f0im   2.5f0 + 6.0621777f0im  -0.5f0 - 7.7942286f0im   0.5f0 - 7.7942286f0im
  -6.5f0 - 4.3301272f0im   6.5f0 - 4.3301272f0im  -7.0f0 - 3.4641016f0im   7.0f0 - 3.4641016f0im
  -7.0f0 + 3.4641016f0im   7.0f0 + 3.4641016f0im  -6.5f0 + 4.3301272f0im   6.5f0 + 4.3301272f0im
Eisenstein primes
Eisenstein primes