eisenstein_integers.sf

#!/usr/bin/ruby

# Simple implementation of Eisenstein integers.

# See also:
#   https://en.wikipedia.org/wiki/Eisenstein_integer

class Eisenstein(a, b) {   # represents: a + b*w, where w = (-1 + i*sqrt(3))/2

    method to_s {
        "Eisenstein(#{a}, #{b})"
    }

    method ==(Eisenstein c) {
        (a == c.a) && (b == c.b)
    }

    method add (Number z) {
        Eisenstein(a+z, b)
    }

    method add (Eisenstein z) {
        var (c,d) = (z.a, z.b)
        Eisenstein(a+c, b+d)
    }

    __CLASS__.alias_method(:add, '+')

    method sub (Eisenstein z) {
        var (c,d) = (z.a, z.b)
        Eisenstein(a-c, b-d)
    }

    __CLASS__.alias_method(:sub, '-')

    method norm {
        a*a - a*b + b*b
    }

    method mul (Number z) {
        Eisenstein(a*z, b*z)
    }

    method mul (Eisenstein z) {
        var (c,d) = (z.a, z.b)
        Eisenstein(a*c - b*d, b*c + a*d - b*d)
    }

    __CLASS__.alias_method(:mul, '*')

    method mod (Number m) {
        Eisenstein(a % m, b % m)
    }

    __CLASS__.alias_method(:mod, '%')

    method pow(Number n) {
        var x = self
        var c = Eisenstein(1, 0)

        for bit in (n.digits(2)) {
            c *= x if bit
            x *= x
        }

        return c
    }

    __CLASS__.alias_method(:pow, '**')

    method powmod(Number n, Number m) {

        var x = self
        var c = Eisenstein(1, 0)

        for bit in (n.digits(2)) {
            (c *= x) %= m if bit        #=
            (x *= x) %= m               #=
        }

        return c
    }
}

# Integers (a,b) such that a^2 - a*b + b^2 give the powers of 7.
with (Eisenstein(2,3)) {|q|
    say 15.of { q.pow(_).a }        #=> [1, 2, -5, -19, 16, 149, 37, -1006, -1265, 5777, 14632, -25807, -128231, 52418, 950035]
    say 15.of { q.pow(_).b }        #=> [0, 3, 3, -18, -39, 87, 360, -249, -2769, -1026, 18357, 25539, -102960, -281733, 438987]
}

var n = (274177-1)
var m = (2**64 + 1)

with (Eisenstein(2,3)) {|q|
    var r = q.powmod(n, m)
    say gcd(r.a-1, m)       #=> 2741177
    say gcd(r.b, m)         #=> 2741177
}