Jump to content

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 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ω|² = 
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


ALGOL 68

Translation of: FreeBASIC
Translation of: Wren

Algol 68 has a standard Complex mode (type) called COMPL - this is used here. This sample doesn't produce the stretch goal plot. The output agrees with the Wren sample.

BEGIN # Eisenstein primes - translated from the FreeBASIC and Wren samples   #

    PROC is prime = ( INT n )BOOL:
         BEGIN
            BOOL result := n = 2 OR ( ODD n AND n > 1 );
            FOR k FROM 3 BY 2 WHILE result AND k * k <= n DO result := n MOD k /= 0 OD;
            result
         END # is prime # ;

    MODE  EISENSTEIN = STRUCT( INT a, b, norm, COMPL n );
    REAL  root3      = sqrt( 3 );
    REAL  half root3 = root3 / 2;

    COMPL omega      = ( -0.5, half root3 );

    PRIO NEWEISENSTEIN = 1, =:= = 1;
    OP   NEWEISENSTEIN = ( INT a, b )EISENSTEIN:
         ( a
         , b
         , ( a * a ) - ( a * b ) + ( b * b )
         , a + omega * b
         );

    OP   ISEISENSTEINPRIME = ( EISENSTEIN e )BOOL:
         IF INT a = a OF e, b = b OF e;
            a = 0 OR b = 0 OR a = b 
         THEN
             INT c = IF ABS a > ABS b THEN ABS a ELSE ABS b FI;
             is prime( c ) AND ( c MOD 3 ) = 2
         ELSE
             is prime( norm OF e )
         FI;
    OP   < = ( EISENSTEIN a, b )BOOL:
         (      INT norm a = norm OF a, norm b = norm OF b;
                norm a < norm b
         OR (   norm a = norm b AND im OF n OF a < im OF n OF b )
         OR (   norm a = norm b AND im OF n OF a = im OF n OF b AND re OF n OF a < re OF n OF b )
         );
    OP   > = ( EISENSTEIN a, b )BOOL:
         (      INT norm a = norm OF a, norm b = norm OF b;
                norm a > norm b
         OR (   norm a = norm b AND im OF n OF a > im OF n OF b )
         OR (   norm a = norm b AND im OF n OF a = im OF n OF b AND re OF n OF a > re OF n OF b )
         );
    OP   =:= = ( REF EISENSTEIN a, b )VOID: BEGIN EISENSTEIN t = a; a := b; b := t END;

    PROC sort eisenstein = ( REF[]EISENSTEIN primes, INT first, last )VOID:
         BEGIN
             INT        i    := first, j := last;
             EISENSTEIN pivot = primes[ ( first + last ) OVER 2 ];
             WHILE WHILE primes[ i ] < pivot DO i +:= 1 OD;
                   WHILE primes[ j ] > pivot DO j -:= 1 OD;
                   IF i <= j THEN
                       primes[ i ] =:= primes[ j ];
                       i           +:= 1;
                       j           -:= 1
                   FI;
                   i <= j
             DO SKIP OD;
             IF first < j THEN sort eisenstein( primes, first, j    ) FI;
             IF i < last  THEN sort eisenstein( primes, i,     last ) FI
         END # sort eisenstein # ;

    BEGIN
        OP   TOSTRING = ( COMPL n )STRING:
             fixed(     re OF n, -7, 4 ) + IF im OF n < 0 THEN " -" ELSE " +" FI
           + fixed( ABS im OF n, -7, 4 ) + "i"
           ;

        INT max range    = 100;
        INT prime count := 0;
        [ 1 : max range * max range * 4 ]EISENSTEIN primes;
        FOR a FROM - max range TO max range DO
            FOR b FROM - max range TO max range DO
                EISENSTEIN e = a NEWEISENSTEIN b;
                IF ISEISENSTEINPRIME e THEN
                    primes[ prime count +:= 1 ] := e
                FI
            OD
        OD;
        sort eisenstein ( primes, 1, prime count );
        print( ( "First 100 Eisenstein primes nearest zero:", newline ) );
        FOR i TO max range DO
            print( ( TOSTRING n OF primes[ i ], "  " ) );
            IF i MOD 4 = 0 THEN print( ( newline ) ) FI
        OD

    END
END
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

FreeBASIC

'#include "isprime.bas"

Type Complex
    real As Double
    imag As Double
End Type

Type Eisenstein
    a As Integer
    b As Integer
    norm As Integer
    n As Complex
End Type

Dim Shared As Complex OMEGA = Type(-0.5, Sqr(3) * 0.5)

' Eisenstein functions
Function EisensteinNew(a As Integer, b As Integer) As Eisenstein
    Dim As Eisenstein e
    e.a = a
    e.b = b
    e.norm = (a * a) - (a * b) + (b * b)
    e.n.real = a - 0.5 * b
    e.n.imag = b * Sqr(3) * 0.5
    Return e
End Function

Function EisensteinIsPrime(e As Eisenstein) As Integer
    If e.a = 0 Or e.b = 0 Or e.a = e.b Then
        Dim As Integer c = Iif(Abs(e.a) > Abs(e.b), Abs(e.a), Abs(e.b))
        Return IsPrime(c) Andalso (c Mod 3) = 2
    End If
    Return IsPrime(e.norm)
End Function

Sub SortEisenstein(primes() As Eisenstein, first As Integer, last As Integer)
    Dim As Integer i = first, j = last
    Dim As Eisenstein pivot = primes((first + last) \ 2)
    
    Do
        While (primes(i).norm < pivot.norm) Or _
            (primes(i).norm = pivot.norm And primes(i).n.imag < pivot.n.imag) Or _
            (primes(i).norm = pivot.norm And primes(i).n.imag = pivot.n.imag And primes(i).n.real < pivot.n.real)
            i += 1
        Wend
        
        While (primes(j).norm > pivot.norm) Or _
            (primes(j).norm = pivot.norm And primes(j).n.imag > pivot.n.imag) Or _
            (primes(j).norm = pivot.norm And primes(j).n.imag = pivot.n.imag And primes(j).n.real > pivot.n.real)
            j -= 1
        Wend
        
        If i <= j Then
            Swap primes(i), primes(j)
            i += 1
            j -= 1
        End If
    Loop Until i > j
    
    If first < j Then SortEisenstein(primes(), first, j)
    If i < last Then SortEisenstein(primes(), i, last)
End Sub

' Main program
Dim As Integer maxRange = 100
Dim As Integer primeCount = 0
Dim As Eisenstein primes(maxRange * maxRange * 4)

For a As Integer = -maxRange To maxRange
    For b As Integer = -maxRange To maxRange
        Dim As Eisenstein e = EisensteinNew(a, b)
        If EisensteinIsPrime(e) Then
            primes(primeCount) = e
            primeCount += 1
        End If
    Next
Next

' Try to replicate Julia sort order for easy comparison
SortEisenstein(primes(), 0, primeCount - 1)

' Display first 100 to terminal
Print "First 100 Eisenstein primes nearest zero:"
For i As Integer = 0 To 99
    Print Using " ##.#### + ##.####i "; primes(i).n.real; primes(i).n.imag;
    If (i + 1) Mod 4 = 0 Then Print
Next

' Generate points array for plotting
Type Point2D
    x As Double
    y As Double
End Type

Function GeneratePlotPoints(primes() As Eisenstein, count As Integer) As Point2D Ptr
    Dim As Point2D Ptr points = Callocate(count * Sizeof(Point2D))
    
    For i As Integer = 0 To count - 1
        points[i].x = primes(i).n.real
        points[i].y = primes(i).n.imag
    Next
    
    Return points
End Function

Dim As Point2D Ptr plotPoints = GeneratePlotPoints(primes(), primeCount)

' Screen setup for plotting
Screenres 1000, 600, 32
Windowtitle "Eisenstein primes with norm <= 100"
Dim As Double scaleX = 3
Dim As Double scaleY = 3
Dim As Double offsetX = 520
Dim As Double offsetY = 280

' Draw axes
Line (offsetX -480, 20)-(offsetX -480, offsetY +280), Rgb(128,128,128) ' Y axis
Line (offsetX -480, offsetY +280)-(980, offsetY +280), Rgb(128,128,128) ' X axis

' Draw scale values
For i As Integer = -50 To 200 Step 25
    Dim As Integer posY = (offsetY +225) - (i * scaleY)
    If (i-75) > -100 Then
        Draw String (offsetX -510, posY -4), Str(i-75)
        Line (offsetX -482, posY)-(offsetX -478, posY), Rgb(128,128,128)
    End If
Next

For i As Integer = -150 To 300 Step 50
    Dim As Integer posX = (offsetX -450) + (i * scaleX)
    Draw String (posX -10, offsetY +290), Str(i-150)
    Line (posX, offsetY +278)-(posX, offsetY +282), Rgb(128,128,128)
Next

' Plot points with screen coordinate conversion
For i As Integer = 0 To primeCount - 1
    Dim As Integer screenX = offsetX + (primes(i).n.real * scaleX)
    Dim As Integer screenY = offsetY - (primes(i).n.imag * scaleY)
    Circle (screenX, screenY), .1, Rgb(255,0,0)
Next

Sleep
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  -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  -0.5000 +  7.7942i   0.5000 +  7.7942i
 -2.5000 + -7.7942i   2.5000 + -7.7942i  -5.5000 + -6.0622i   5.5000 + -6.0622i
Eisenstein primes
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

jq

Works with jq, the C implementation of jq

Works with gojq, the Go implementation of jq

Adapted from Wren

In this entry, complex numbers are represented as arrays of pairs: [real, complex], as in the jq section on the Complex page. The two functions for adding and multiplying complex numbers presented there are reproduced below so that the program presented here is self-contained.

For the "stretch" task, we assume the availability of a tool such as gnuplot; using gnuplot, a suitable sequence of commands to plot the points produced by `graph` as defined below would be as follows:

reset
set terminal pngcairo
set output "eisenstein-primes.png"
### Complex numbers
def plus(x; y):
    if (x|type) == "number" then
       if  (y|type) == "number" then [ x+y, 0 ]
       else [ x + y[0], y[1]]
       end
    elif (y|type) == "number" then plus(y;x)
    else [ x[0] + y[0], x[1] + y[1] ]
    end;
 
def multiply(x; y):
    if (x|type) == "number" then
       if  (y|type) == "number" then [ x*y, 0 ]
       else [x * y[0], x * y[1]]
       end
    elif (y|type) == "number" then multiply(y;x)
    else [ x[0] * y[0] - x[1] * y[1],  x[0] * y[1] + x[1] * y[0]]
    end;

### Generic utilities

def lpad($len): tostring | ($len - length) as $l | (" " * $l) + .;

# Require $n > 0
def nwise($n):
  def _n: if length <= $n then . else .[:$n] , (.[$n:] | _n) end;
  if $n <= 0 then "nwise: argument should be non-negative" else _n end;

def is_prime:
  . as $n
  | if ($n < 2)         then false
    elif ($n % 2 == 0)  then $n == 2
    elif ($n % 3 == 0)  then $n == 3
    elif ($n % 5 == 0)  then $n == 5
    elif ($n % 7 == 0)  then $n == 7
    elif ($n % 11 == 0) then $n == 11
    elif ($n % 13 == 0) then $n == 13
    elif ($n % 17 == 0) then $n == 17
    elif ($n % 19 == 0) then $n == 19
    else sqrt as $s
    | 23
    | until( . > $s or ($n % . == 0); . + 2)
    | . > $s
    end;

def OMEGA: [-0.5, (3|sqrt * 0.5)];

### Eisenstein numbers and Eisenstein primes

def Eisenstein($a; $b):
   {$a, $b, n: plus( multiply(OMEGA;$b); $a) };

def realEisenstein: .n[0];
def imagEisenstein: .n[1];
def normEisenstein:
   .a *.a - .a * .b + .b * .b ;

# Replicate the Julia sort order for easy comparison
def sortEisenstein:
   sort_by( [ normEisenstein, imagEisenstein, realEisenstein] );

def isPrimeEisenstein:
    if .a == 0 or .b == 0 or .a == .b
    # length ~ abs
    then ([.a, .b] | map(length) | max) as $c  
    | ($c | is_prime) and $c % 3 == 2
    else normEisenstein | is_prime
    end;

# Eisenstein($i;$j) primes for $i and $j in -$n .. $n inclusive
def eprimes($n):
  reduce range (-$n; $n+1) as $a ([];
    reduce range ( -$n; $n+1) as $b (.;
      Eisenstein($a; $b) as $e
      | if $e | isPrimeEisenstein
        then . + [$e]
        else .
        end ));

### The tasks

# pretty-print a complex number
def pp:
  def r: 100 * . | trunc / 100;
  .[2] = (if .[1] < 0 then "-" else "+" end)
  | .[1] |= (if . < 0 then -. else . end)
  | "\(.[0]|r|lpad(5)) \(.[2]) \(.[1]|r|lpad(5))i";

# Display the input array of complex numbers as a table with $n columns
# proceeding row-wise and using pp/0
def row_wise($n):
  nwise($n) | map( pp ) | join("  ");

def listing:
  {eprimes: (eprimes(10) | sortEisenstein) }
  # convert to Complex numbers for easy display
  | .eprimes |= map( .n )
  | "First 100 Eisenstein primes nearest zero:",
    (.eprimes[:100] | row_wise(4) );

def graph:
  eprimes(100)
  | sortEisenstein
  | .[:2000][]
  | .n
  | "\(real(.)) \(imag(.))";

# For a listing of the first 100 Eisenstein primes nearest 0:
listing

# To produce the points for gnuplot:
# graph
Output:

The results of `listing` are shown below. For the graph of the output produced by `graph`, see the graph shown above at J.

First 100 Eisenstein primes nearest zero:
    0 -  1.73i   -1.5 -  0.86i    1.5 -  0.86i   -1.5 +  0.86i
  1.5 +  0.86i      0 +  1.73i     -1 -  1.73i      1 -  1.73i
   -2 +     0i      2 +     0i     -1 +  1.73i      1 +  1.73i
 -0.5 -  2.59i    0.5 -  2.59i     -2 -  1.73i      2 -  1.73i
 -2.5 -  0.86i    2.5 -  0.86i   -2.5 +  0.86i    2.5 +  0.86i
   -2 +  1.73i      2 +  1.73i   -0.5 +  2.59i    0.5 +  2.59i
   -1 -  3.46i      1 -  3.46i   -2.5 -  2.59i    2.5 -  2.59i
 -3.5 -  0.86i    3.5 -  0.86i   -3.5 +  0.86i    3.5 +  0.86i
 -2.5 +  2.59i    2.5 +  2.59i     -1 +  3.46i      1 +  3.46i
 -0.5 -  4.33i    0.5 -  4.33i   -3.5 -  2.59i    3.5 -  2.59i
   -4 -  1.73i      4 -  1.73i     -4 +  1.73i      4 +  1.73i
 -3.5 +  2.59i    3.5 +  2.59i   -0.5 +  4.33i    0.5 +  4.33i
 -2.5 -  4.33i    2.5 -  4.33i     -5 +     0i      5 +     0i
 -2.5 +  4.33i    2.5 +  4.33i     -2 -  5.19i      2 -  5.19i
 -3.5 -  4.33i    3.5 -  4.33i   -5.5 -  0.86i    5.5 -  0.86i
 -5.5 +  0.86i    5.5 +  0.86i   -3.5 +  4.33i    3.5 +  4.33i
   -2 +  5.19i      2 +  5.19i   -0.5 -  6.06i    0.5 -  6.06i
   -5 -  3.46i      5 -  3.46i   -5.5 -  2.59i    5.5 -  2.59i
 -5.5 +  2.59i    5.5 +  2.59i     -5 +  3.46i      5 +  3.46i
 -0.5 +  6.06i    0.5 +  6.06i   -2.5 -  6.06i    2.5 -  6.06i
   -4 -  5.19i      4 -  5.19i   -6.5 -  0.86i    6.5 -  0.86i
 -6.5 +  0.86i    6.5 +  0.86i     -4 +  5.19i      4 +  5.19i
 -2.5 +  6.06i    2.5 +  6.06i   -0.5 -  7.79i    0.5 -  7.79i
 -6.5 -  4.33i    6.5 -  4.33i     -7 -  3.46i      7 -  3.46i
   -7 +  3.46i      7 +  3.46i   -6.5 +  4.33i    6.5 +  4.33i

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
Eisenstein primes

Nim

Translation of: Wren
import std/[algorithm, complex, math, strformat]

import gnuplot

func isPrime(n: Natural): bool =
  if n < 2: return false
  if (n and 1) == 0: return n == 2
  if n mod 3 == 0: return n == 3
  var k = 5
  var delta = 2
  while k * k <= n:
    if n mod k == 0: return false
    inc k, delta
    delta = 6 - delta
  result = true


### Eisenstein definition.

const ω = complex(-0.5, sqrt(3.0) * 0.5)

type Eisenstein = object
  a: int
  b: int
  n: Complex64

func initEisenstein(a, b: int): Eisenstein =
  ## Initialize an Eisenstein number.
  Eisenstein(a: a, b: b, n: a.toFloat + b.toFloat * ω)

template re(e: Eisenstein): float = e.n.re
template im(e: Eisenstein): float = e.n.im

func norm(e: Eisenstein): int =
  ## return the norm of an Eisenstein number.
  e.a * e.a - e.a * e.b + e.b * e.b

func isPrime(e: Eisenstein): bool =
  ## Return true if an Eisenstein number is prime.
  if e.a == 0 or e.b == 0 or e.a == e.b:
    let c = max(abs(e.a), abs(e.b))
    result = c.isPrime and c mod 3 == 2
  else:
    result = e.norm.isPrime

func `$`(e: Eisenstein): string =
  ## Return a string representation of an Eisenstein number.
  let (sign, im) = if e.im >= 0: ('+', e.im) else: ('-', -e.im)
  result =  &"{e.re:7.4f} {sign} {im:6.4f}i"


### Find Eisenstein primes.

var eprimes: seq[Eisenstein]
for a in -100..100:
  for b in -100..100:
    let e = initEisenstein(a, b)
    if e.isPrime: eprimes.add e

# Try to replicate Wren sort order for easy comparison.
eprimes.sort(proc (e1, e2: Eisenstein): int =
               result = cmp(e1.norm, e2.norm)
               if result == 0:
                 result = cmp(e1.im, e2.im)
                 if result == 0:
                   result = cmp(e1.re, e2.re)
            )

# Display first 100 Eisenstein primes to terminal.
echo "First 100 Eisenstein primes nearest zero:"
for i in 0..99:
  stdout.write eprimes[i]
  stdout.write if i mod 4 == 3: "\n" else: "  "

# Generate points for the plot.
var x, y: seq[float]
for e in eprimes:
  x.add e.re
  y.add e.im

withGnuPlot:
  cmd "set size ratio -1"
  plot(x, y, "Eisenstein primes", "with dots lw 2")
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
Eisenstein primes
Eisenstein primes

Perl

Translation of: Raku
use v5.36;
use Math::AnyNum <pi mod max complex reals is_prime>;

my $omega = exp ( complex(0,2) * pi/3 ); my @E;

sub    norm (@p)     { $p[0]**2 - $p[0]*$p[1] + $p[1]**2 }
sub display (@p)     { sprintf '%+8.4f%+8.4fi', reals($p[0] + $omega*$p[1]) }
sub       X ($a, $b) { my @p; for my $x ($a..$b) { for my $y ($a..$b) { push @p, [$x, $y] } } @p }
sub   table ($c, @V) { my $t = $c * (my $w = 1 + max map { length } @V); ( sprintf( ('%'.$w.'s')x@V, @V) ) =~ s/.{1,$t}\K/\n/gr }

for (X -10, 10) {
    my($a,$b) = @$_;
    my $c = max abs($a), abs($b);
    push @E, [@$_] if ((0==$a or 0==$b or $a==$b) and is_prime $c and 2 == mod $c,3) or is_prime norm @$_
}
say table 4, (map { display @$_ } sort { norm(@$a) <=> norm(@$b) } @E)[0..99];
Output:
  -1.5000 -0.8660i  -0.0000 -1.7321i  -1.5000 +0.8660i  +1.5000 -0.8660i
  +0.0000 +1.7321i  +1.5000 +0.8660i  -1.0000 -1.7321i  -2.0000 +0.0000i
  +1.0000 -1.7321i  -1.0000 +1.7321i  +2.0000 +0.0000i  +1.0000 +1.7321i
  -2.0000 -1.7321i  -2.5000 -0.8660i  -0.5000 -2.5981i  -2.5000 +0.8660i
  +0.5000 -2.5981i  -2.0000 +1.7321i  +2.0000 -1.7321i  -0.5000 +2.5981i
  +2.5000 -0.8660i  +0.5000 +2.5981i  +2.5000 +0.8660i  +2.0000 +1.7321i
  -2.5000 -2.5981i  -3.5000 -0.8660i  -1.0000 -3.4641i  -3.5000 +0.8660i
  +1.0000 -3.4641i  -2.5000 +2.5981i  +2.5000 -2.5981i  -1.0000 +3.4641i
  +3.5000 -0.8660i  +1.0000 +3.4641i  +3.5000 +0.8660i  +2.5000 +2.5981i
  -3.5000 -2.5981i  -4.0000 -1.7321i  -0.5000 -4.3301i  -4.0000 +1.7321i
  +0.5000 -4.3301i  -3.5000 +2.5981i  +3.5000 -2.5981i  -0.5000 +4.3301i
  +4.0000 -1.7321i  +0.5000 +4.3301i  +4.0000 +1.7321i  +3.5000 +2.5981i
  -2.5000 -4.3301i  -5.0000 +0.0000i  +2.5000 -4.3301i  -2.5000 +4.3301i
  +5.0000 +0.0000i  +2.5000 +4.3301i  -3.5000 -4.3301i  -5.5000 -0.8660i
  -2.0000 -5.1962i  -5.5000 +0.8660i  +2.0000 -5.1962i  -3.5000 +4.3301i
  +3.5000 -4.3301i  -2.0000 +5.1962i  +5.5000 -0.8660i  +2.0000 +5.1962i
  +5.5000 +0.8660i  +3.5000 +4.3301i  -5.0000 -3.4641i  -5.5000 -2.5981i
  -0.5000 -6.0622i  -5.5000 +2.5981i  +0.5000 -6.0622i  -5.0000 +3.4641i
  +5.0000 -3.4641i  -0.5000 +6.0622i  +5.5000 -2.5981i  +0.5000 +6.0622i
  +5.5000 +2.5981i  +5.0000 +3.4641i  -4.0000 -5.1962i  -6.5000 -0.8660i
  -2.5000 -6.0622i  -6.5000 +0.8660i  +2.5000 -6.0622i  -4.0000 +5.1962i
  +4.0000 -5.1962i  -2.5000 +6.0622i  +6.5000 -0.8660i  +2.5000 +6.0622i
  +6.5000 +0.8660i  +4.0000 +5.1962i  -6.5000 -4.3301i  -7.0000 -3.4641i
  -0.5000 -7.7942i  -7.0000 +3.4641i  +0.5000 -7.7942i  -6.5000 +4.3301i
  +6.5000 -4.3301i  -0.5000 +7.7942i  +7.0000 -3.4641i  +0.5000 +7.7942i

Phix

Library: Phix/xpGUI
requires("1.0.3")
include complex.e
constant OMEGA = {-0.5, sqrt(3)*0.5}

// try to replicate Wren sort order for easy comparison
enum NORM, IMAG, REAL, IS_PRIME, ELEN=$
function new_Eisenstein(integer a,b)
    atom {real,imag} = complex_add(complex_mul(OMEGA,b),a)
    integer norm = a*a-a*b+b*b, c = max(abs(a),abs(b))
    bool p = iff(a=0 or b=0 or a=b?is_prime(c) and remainder(c,3)=2
                                  :is_prime(norm))
    return {norm,imag,real,p} -- nb in [NORM..IS_PRIME] order
end function

function Eisenstein()
    sequence eprimes = {}
    for a=-100 to 100 do
        for b=-100 to 100 do
            sequence e = new_Eisenstein(a, b)
            if e[IS_PRIME] then
                eprimes = append(eprimes,e)
            end if
        end for
    end for
    eprimes = sort(eprimes)
    sequence real = vslice(eprimes,REAL),
             imag = vslice(eprimes,IMAG),
             f100 = repeat("",100)
    for i=1 to 100 do -- convert for display
        integer pm = iff(imag[i]>=0?'+':'-')
        f100[i] = sprintf("%7.4f %c%7.4fi",{real[i],pm,abs(imag[i])})
    end for
    return {f100,real,imag}
end function

sequence {f100,real,imag} = Eisenstein()
printf(1,"First 100 Eisenstein primes nearest zero:\n%s\n",join_by(f100,1,4,"  "))

include xpGUI.e
function get_data(gdx graph) return {{real,imag}} end function
gdx graph = gGraph(get_data,"XMIN=-150,XMAX=150,YMIN=-100,YMAX=100"),
    dlg = gDialog(graph,"Eisenstein primes","SIZE=392x290")
gSetAttribute(graph,"GTITLE","with norm <= 100  (%d points)",{length(real)})
gSetAttributes(graph,"XTICK=50,YTICK=25,MARKSTYLE=DOT,GRID=NO")
gShow(dlg)
gMainLoop()

Output same as Wren.

Raku

my \ω = exp 2i × π/3;

sub norm    (@p) { @p[0]² - @p[0@p[1] + @p[1]² }
sub display (@p) { (@p[0] + ω×@p[1]).reals».fmt('%+8.4f').join ~ 'i' }

my @E = gather (-10..10 X -10..10).map: -> (\a,\b) {
    take (a,b) if 0 == a|b || a == b ?? (.is-prime and 2 == $_ mod 3 given (a,b)».abs.max) !! norm((a,b)).is-prime
}

(@E.sort: *.&norm).head(100).map(*.&display).batch(4).join("\n").say;
Output:
 -1.5000 -0.8660i  -0.0000 -1.7321i  -1.5000 +0.8660i  +1.5000 -0.8660i
 +0.0000 +1.7321i  +1.5000 +0.8660i  -1.0000 -1.7321i  -2.0000 +0.0000i
 +1.0000 -1.7321i  -1.0000 +1.7321i  +2.0000 +0.0000i  +1.0000 +1.7321i
 -2.0000 -1.7321i  -2.5000 -0.8660i  -0.5000 -2.5981i  -2.5000 +0.8660i
 +0.5000 -2.5981i  -2.0000 +1.7321i  +2.0000 -1.7321i  -0.5000 +2.5981i
 +2.5000 -0.8660i  +0.5000 +2.5981i  +2.5000 +0.8660i  +2.0000 +1.7321i
 -2.5000 -2.5981i  -3.5000 -0.8660i  -1.0000 -3.4641i  -3.5000 +0.8660i
 +1.0000 -3.4641i  -2.5000 +2.5981i  +2.5000 -2.5981i  -1.0000 +3.4641i
 +3.5000 -0.8660i  +1.0000 +3.4641i  +3.5000 +0.8660i  +2.5000 +2.5981i
 -3.5000 -2.5981i  -4.0000 -1.7321i  -0.5000 -4.3301i  -4.0000 +1.7321i
 +0.5000 -4.3301i  -3.5000 +2.5981i  +3.5000 -2.5981i  -0.5000 +4.3301i
 +4.0000 -1.7321i  +0.5000 +4.3301i  +4.0000 +1.7321i  +3.5000 +2.5981i
 -2.5000 -4.3301i  -5.0000 +0.0000i  +2.5000 -4.3301i  -2.5000 +4.3301i
 +5.0000 +0.0000i  +2.5000 +4.3301i  -3.5000 -4.3301i  -5.5000 -0.8660i
 -2.0000 -5.1962i  -5.5000 +0.8660i  +2.0000 -5.1962i  -3.5000 +4.3301i
 +3.5000 -4.3301i  -2.0000 +5.1962i  +5.5000 -0.8660i  +2.0000 +5.1962i
 +5.5000 +0.8660i  +3.5000 +4.3301i  -5.0000 -3.4641i  -5.5000 -2.5981i
 -0.5000 -6.0622i  -5.5000 +2.5981i  +0.5000 -6.0622i  -5.0000 +3.4641i
 +5.0000 -3.4641i  -0.5000 +6.0622i  +5.5000 -2.5981i  +0.5000 +6.0622i
 +5.5000 +2.5981i  +5.0000 +3.4641i  -4.0000 -5.1962i  -6.5000 -0.8660i
 -2.5000 -6.0622i  -6.5000 +0.8660i  +2.5000 -6.0622i  -4.0000 +5.1962i
 +4.0000 -5.1962i  -2.5000 +6.0622i  +6.5000 -0.8660i  +2.5000 +6.0622i
 +6.5000 +0.8660i  +4.0000 +5.1962i  -6.5000 -4.3301i  -7.0000 -3.4641i
 -0.5000 -7.7942i  -7.0000 +3.4641i  +0.5000 -7.7942i  -6.5000 +4.3301i
 +6.5000 -4.3301i  -0.5000 +7.7942i  +7.0000 -3.4641i  +0.5000 +7.7942i

Sidef

Translation of: Raku
class Eisenstein(a, b, w = (-1 + sqrt(3).i)/2) {
    method norm {
        a**2 - a*b + b**2
    }

    method to_s {
        sprintf('%+8.4f%+8.4fi', reals(a + b*w))
    }
}

var E = []

for e in (-10..10 ~X -10..10 -> map_2d {|x,y| Eisenstein(x,y) }) {
    var c = [e.a,e.b].map{.abs}.max
    if (
        ((0 ~~ [e.a, e.b]) || (e.a == e.b)) ?
        (c.is_congruent(2,3) && c.is_prime) : e.norm.is_prime
    ) {
        E << e
    }
}

E.sort_by { .norm }.first(100).slices(4).each {|s|
    say s.join('  ')
}
Output:
 -1.5000 -0.8660i   +0.0000 -1.7321i   -1.5000 +0.8660i   +1.5000 -0.8660i
 +0.0000 +1.7321i   +1.5000 +0.8660i   -1.0000 -1.7321i   -2.0000 +0.0000i
 +1.0000 -1.7321i   -1.0000 +1.7321i   +2.0000 +0.0000i   +1.0000 +1.7321i
 -2.0000 -1.7321i   -2.5000 -0.8660i   -0.5000 -2.5981i   -2.5000 +0.8660i
 +0.5000 -2.5981i   -2.0000 +1.7321i   +2.0000 -1.7321i   -0.5000 +2.5981i
 +2.5000 -0.8660i   +0.5000 +2.5981i   +2.5000 +0.8660i   +2.0000 +1.7321i
 -2.5000 -2.5981i   -3.5000 -0.8660i   -1.0000 -3.4641i   -3.5000 +0.8660i
 +1.0000 -3.4641i   -2.5000 +2.5981i   +2.5000 -2.5981i   -1.0000 +3.4641i
 +3.5000 -0.8660i   +1.0000 +3.4641i   +3.5000 +0.8660i   +2.5000 +2.5981i
 -3.5000 -2.5981i   -4.0000 -1.7321i   -0.5000 -4.3301i   -4.0000 +1.7321i
 +0.5000 -4.3301i   -3.5000 +2.5981i   +3.5000 -2.5981i   -0.5000 +4.3301i
 +4.0000 -1.7321i   +0.5000 +4.3301i   +4.0000 +1.7321i   +3.5000 +2.5981i
 -2.5000 -4.3301i   -5.0000 +0.0000i   +2.5000 -4.3301i   -2.5000 +4.3301i
 +5.0000 +0.0000i   +2.5000 +4.3301i   -3.5000 -4.3301i   -5.5000 -0.8660i
 -2.0000 -5.1962i   -5.5000 +0.8660i   +2.0000 -5.1962i   -3.5000 +4.3301i
 +3.5000 -4.3301i   -2.0000 +5.1962i   +5.5000 -0.8660i   +2.0000 +5.1962i
 +5.5000 +0.8660i   +3.5000 +4.3301i   -5.0000 -3.4641i   -5.5000 -2.5981i
 -0.5000 -6.0622i   -5.5000 +2.5981i   +0.5000 -6.0622i   -5.0000 +3.4641i
 +5.0000 -3.4641i   -0.5000 +6.0622i   +5.5000 -2.5981i   +0.5000 +6.0622i
 +5.5000 +2.5981i   +5.0000 +3.4641i   -4.0000 -5.1962i   -6.5000 -0.8660i
 -2.5000 -6.0622i   -6.5000 +0.8660i   +2.5000 -6.0622i   -4.0000 +5.1962i
 +4.0000 -5.1962i   -2.5000 +6.0622i   +6.5000 -0.8660i   +2.5000 +6.0622i
 +6.5000 +0.8660i   +4.0000 +5.1962i   -6.5000 -4.3301i   -7.0000 -3.4641i
 -0.5000 -7.7942i   -7.0000 +3.4641i   +0.5000 -7.7942i   -6.5000 +4.3301i
 +6.5000 -4.3301i   -0.5000 +7.7942i   +7.0000 -3.4641i   +0.5000 +7.7942i

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 
Cookies help us deliver our services. By using our services, you agree to our use of cookies.