Eisenstein primes

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: ${\displaystyle \omega ={\frac {-1+i{\sqrt {3}}}{2}}=e^{i2\pi /3}}$

As with 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ω|² = $a^{2}-ab+b^{2}$

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

• [Wikipedia entry for Eisenstein_integer]
• [Observable entry on Eisenstein primes]

J

Implementation:

eisensteinprimes=: {{
  rY=. >.1.5%:y
  p1=. ,(w^i.3)*/(#~ 2= 3|]) p:i.rY
  'a b'=. |:(2#rY)#:I.,1 p: {{(x*x)+(y*y)-x*y}}"0/~i.rY
  y{.(/: *:@|)p1,(,-)(a+b*w),a+b*-w
}}

Task example (and stretch - taking the stretch goal in a minimalist literal fashion):

   20 5$eisensteinprimes 100
     0j1.73205   1.5j0.866025     0j_1.73205 _1.5j_0.866025   _1j_1.73205
             2     _1j1.73205   _0.5j2.59808    0.5j2.59808  2.5j0.866025
     2j1.73205     2j_1.73205  2.5j_0.866025   0.5j_2.59808 _0.5j_2.59808
_2.5j_0.866025    _2j_1.73205     _2j1.73205  _2.5j0.866025     _1j3.4641
      1j3.4641      1j_3.4641     _1j_3.4641   3.5j0.866025   2.5j2.59808
  2.5j_2.59808  3.5j_0.866025 _3.5j_0.866025  _2.5j_2.59808  _2.5j2.59808
 _3.5j0.866025   _0.5j4.33013    0.5j4.33013    3.5j2.59808  3.5j_2.59808
  0.5j_4.33013  _0.5j_4.33013  _3.5j_2.59808   _3.5j2.59808     4j1.73205
    4j_1.73205    _4j_1.73205     _4j1.73205      3j_3.4641     _3j3.4641
 4.5j_0.866025  _4.5j0.866025  _2.5j_4.33013              5  _2.5j4.33013
    _2j5.19615      2j5.19615   5.5j0.866025    3.5j4.33013    2j_5.19615
   _2j_5.19615 _5.5j_0.866025  _3.5j_4.33013   _0.5j6.06218   0.5j6.06218
      5j3.4641      5j_3.4641   0.5j_6.06218  _0.5j_6.06218    _5j_3.4641
     _5j3.4641    5.5j2.59808   5.5j_2.59808  _5.5j_2.59808  _5.5j2.59808
  4.5j_4.33013   _4.5j4.33013     6j_1.73205     _6j1.73205  _2.5j6.06218
   2.5j6.06218   6.5j0.866025      4j5.19615     4j_5.19615 6.5j_0.866025
  2.5j_6.06218  _2.5j_6.06218 _6.5j_0.866025    _4j_5.19615    _4j5.19615
 _6.5j0.866025   5.5j_4.33013   6.5j_2.59808   _5.5j4.33013  _6.5j2.59808
  4.5j_6.06218  7.5j_0.866025   _4.5j6.06218  _7.5j0.866025  _0.5j7.79423
   0.5j7.79423    6.5j4.33013   0.5j_7.79423  _0.5j_7.79423 _6.5j_4.33013

   require'plot'
   'marker; markersize 0.3' plot eisensteinprimes 2000

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<:Integer} = new{T}(a, zero(T))
    Eisenstein(a::Integer, b::Integer) = 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!(filter(is_eisenstein_prime, arr), lt = (x, y) -> norm(x) < norm(y)))
    for (i, c) in enumerate(eprimes)
        if i <= printlimit
            print(lpad(round(Complex(c), digits = 4), 18), i % 5 == 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.0 - 1.7321im    -1.5 - 0.866im     1.5 - 0.866im    -1.5 + 0.866im     1.5 + 0.866im
    0.0 + 1.7321im   -1.0 - 1.7321im    1.0 - 1.7321im      -2.0 + 0.0im       2.0 + 0.0im
   -1.0 + 1.7321im    1.0 + 1.7321im   -0.5 - 2.5981im    0.5 - 2.5981im   -2.0 - 1.7321im
    2.0 - 1.7321im    -2.5 - 0.866im     2.5 - 0.866im    -2.5 + 0.866im     2.5 + 0.866im
   -2.0 + 1.7321im    2.0 + 1.7321im   -0.5 + 2.5981im    0.5 + 2.5981im   -1.0 - 3.4641im
    1.0 - 3.4641im   -2.5 - 2.5981im    2.5 - 2.5981im    -3.5 - 0.866im     3.5 - 0.866im
    -3.5 + 0.866im     3.5 + 0.866im   -2.5 + 2.5981im    2.5 + 2.5981im   -1.0 + 3.4641im
    1.0 + 3.4641im   -0.5 - 4.3301im    0.5 - 4.3301im   -3.5 - 2.5981im    3.5 - 2.5981im
   -4.0 - 1.7321im    4.0 - 1.7321im   -4.0 + 1.7321im    4.0 + 1.7321im   -3.5 + 2.5981im
    3.5 + 2.5981im   -0.5 + 4.3301im    0.5 + 4.3301im   -2.5 - 4.3301im    2.5 - 4.3301im
      -5.0 + 0.0im       5.0 + 0.0im   -2.5 + 4.3301im    2.5 + 4.3301im   -2.0 - 5.1962im
    2.0 - 5.1962im   -3.5 - 4.3301im    3.5 - 4.3301im    -5.5 - 0.866im     5.5 - 0.866im
    -5.5 + 0.866im     5.5 + 0.866im   -3.5 + 4.3301im    3.5 + 4.3301im   -2.0 + 5.1962im
    2.0 + 5.1962im   -0.5 - 6.0622im    0.5 - 6.0622im   -5.0 - 3.4641im    5.0 - 3.4641im
   -5.5 - 2.5981im    5.5 - 2.5981im   -5.5 + 2.5981im    5.5 + 2.5981im   -5.0 + 3.4641im
    5.0 + 3.4641im   -0.5 + 6.0622im    0.5 + 6.0622im   -2.5 - 6.0622im    2.5 - 6.0622im
   -4.0 - 5.1962im    4.0 - 5.1962im    -6.5 - 0.866im     6.5 - 0.866im    -6.5 + 0.866im
     6.5 + 0.866im   -4.0 + 5.1962im    4.0 + 5.1962im   -2.5 + 6.0622im    2.5 + 6.0622im
   -0.5 - 7.7942im    0.5 - 7.7942im   -6.5 - 4.3301im    6.5 - 4.3301im   -7.0 - 3.4641im
    7.0 - 3.4641im   -7.0 + 3.4641im    7.0 + 3.4641im   -6.5 + 4.3301im    6.5 + 4.3301im

Eisenstein primes

Wren

Library: DOME

Library: Wren-plot

Library: Wren-iterate

Library: Wren-complex

Library: Wren-math

Library: Wren-fmt

import "dome" for Window
import "graphics" for Canvas, Color
import "./plot" for Axes
import "./iterate" for Stepped
import "./complex" for Complex
import "./math2" for Math, Int
import "./fmt" for Fmt

var OMEGA = Complex.new(-0.5, 3.sqrt * 0.5)

class Eisenstein {
    construct new(a, b) {
        _a = a
        _b = b
        _n = OMEGA * b + a
    }

    a { _a }
    b { _b }
    n { _n }

    real { _n.real }
    imag { _n.imag }
    norm { _a *_a - _a * _b + _b * _b }

    isPrime {
        if (_a == 0 || _b == 0 || _a == _b) {
            var c = Math.max(_a.abs, _b.abs)
            return Int.isPrime(c) && c % 3 == 2
        }
        return Int.isPrime(norm)
    }

    toString { _n.toString }
}

var eprimes = []
for (a in -100..100) {
    for (b in -100..100) {
        var e = Eisenstein.new(a, b)
        if (e.isPrime) eprimes.add(e)
    }
}

// try to replicate Julia sort order for easy comparison
eprimes.sort { |e1, e2|
    if (e1.norm < e2.norm) return true
    if (e1.norm == e2.norm) {
        if (e1.imag < e2.imag) return true
        if (e1.imag == e2.imag) return e1.real < e2.real
        return false
    }
    return false
}

// convert to Complex numbers for easy display
eprimes = eprimes.map { |e| e.n }

// display first 100 to terminal
System.print("First 100 Eisenstein primes nearest zero:")
Fmt.tprint("$ 6.4z ", eprimes.take(100), 4)

// generate points for the plot
var Pts = eprimes.map { |e| [e.real, e.imag] }.toList

class Main {
    construct new() {
        Window.title = "Eisenstein primes with norm <= 100  (%(Pts.count) points)"
        Canvas.resize(1000, 600)
        Window.resize(1000, 600)
        Canvas.cls(Color.white)
        var axes = Axes.new(100, 500, 800, 400, -160..160, -100..100)
        axes.draw(Color.black, 2)
        var xMarks = Stepped.new(-150..150, 50)
        var yMarks = Stepped.new(-75..75, 25)
        axes.mark(xMarks, yMarks, Color.black, 2)
        axes.label(xMarks, yMarks, Color.black, 2, Color.black)
        axes.plot(Pts, Color.black, "·") // uses interpunct character 0xb7
    }

    init() {}

    update() {}

    draw(alpha) {}
}

var Game = Main.new()

Output:

Terminal output:

First 100 Eisenstein primes nearest zero:
 0.0000 - 1.7321i  -1.5000 - 0.8660i   1.5000 - 0.8660i  -1.5000 + 0.8660i 
 1.5000 + 0.8660i   0.0000 + 1.7321i  -1.0000 - 1.7321i   1.0000 - 1.7321i 
-2.0000 + 0.0000i   2.0000 + 0.0000i  -1.0000 + 1.7321i   1.0000 + 1.7321i 
-0.5000 - 2.5981i   0.5000 - 2.5981i  -2.0000 - 1.7321i   2.0000 - 1.7321i 
-2.5000 - 0.8660i   2.5000 - 0.8660i  -2.5000 + 0.8660i   2.5000 + 0.8660i 
-2.0000 + 1.7321i   2.0000 + 1.7321i  -0.5000 + 2.5981i   0.5000 + 2.5981i 
-1.0000 - 3.4641i   1.0000 - 3.4641i  -2.5000 - 2.5981i   2.5000 - 2.5981i 
-3.5000 - 0.8660i   3.5000 - 0.8660i  -3.5000 + 0.8660i   3.5000 + 0.8660i 
-2.5000 + 2.5981i   2.5000 + 2.5981i  -1.0000 + 3.4641i   1.0000 + 3.4641i 
-0.5000 - 4.3301i   0.5000 - 4.3301i  -3.5000 - 2.5981i   3.5000 - 2.5981i 
-4.0000 - 1.7321i   4.0000 - 1.7321i  -4.0000 + 1.7321i   4.0000 + 1.7321i 
-3.5000 + 2.5981i   3.5000 + 2.5981i  -0.5000 + 4.3301i   0.5000 + 4.3301i 
-2.5000 - 4.3301i   2.5000 - 4.3301i  -5.0000 + 0.0000i   5.0000 + 0.0000i 
-2.5000 + 4.3301i   2.5000 + 4.3301i  -2.0000 - 5.1962i   2.0000 - 5.1962i 
-3.5000 - 4.3301i   3.5000 - 4.3301i  -5.5000 - 0.8660i   5.5000 - 0.8660i 
-5.5000 + 0.8660i   5.5000 + 0.8660i  -3.5000 + 4.3301i   3.5000 + 4.3301i 
-2.0000 + 5.1962i   2.0000 + 5.1962i  -0.5000 - 6.0622i   0.5000 - 6.0622i 
-5.0000 - 3.4641i   5.0000 - 3.4641i  -5.5000 - 2.5981i   5.5000 - 2.5981i 
-5.5000 + 2.5981i   5.5000 + 2.5981i  -5.0000 + 3.4641i   5.0000 + 3.4641i 
-0.5000 + 6.0622i   0.5000 + 6.0622i  -2.5000 - 6.0622i   2.5000 - 6.0622i 
-4.0000 - 5.1962i   4.0000 - 5.1962i  -6.5000 - 0.8660i   6.5000 - 0.8660i 
-6.5000 + 0.8660i   6.5000 + 0.8660i  -4.0000 + 5.1962i   4.0000 + 5.1962i 
-2.5000 + 6.0622i   2.5000 + 6.0622i  -0.5000 - 7.7942i   0.5000 - 7.7942i 
-6.5000 - 4.3301i   6.5000 - 4.3301i  -7.0000 - 3.4641i   7.0000 - 3.4641i 
-7.0000 + 3.4641i   7.0000 + 3.4641i  -6.5000 + 4.3301i   6.5000 + 4.3301i