Dual EC Backdoor in Python

Functions.py

#https://leetarxiv.substack.com/p/dual-ec-backdoor-coding-guide
from fastecdsa.curve import P256
from fastecdsa.point import Point
from random import randint

class DualEC():
    def __init__(self, seed, P, Q):
        self.seed = seed # Initial integer state of RNG
        self.P = P # First elliptic curve point (public parameter)
        self.Q = Q # Second elliptic curve point (could be maliciously chosen)

    def GenerateBits(self):
        t = self.seed
        s = (t * self.P).x # Multiply seed by P and take x-coordinate to get new state
        self.seed = s # Update internal state
        r = (s * self.Q).x # Multiply new state by Q and take x-coordinate to get output
        return r & (2**(8 * 30) - 1)  # Truncate to 30 bytes (240 bits)
    

def ModularInverse(a, m):
    # only works if m is prime (due to Euler's Theorem)
    return pow(a, m-2, m)

def GenerateBackdoorConstants():
    P = P256.G  # set P to the curve base
    d = randint(2, P256.q)  # pick a number random number in the field
    e = ModularInverse(d, P256.q)  # find inverse of the number in the field
    Q = e * P  #This is exponentiation in P256
    return P, Q, d

def p256_mod_sqrt(c):
    # only works for field P256 is over
    p = 0xffffffff00000001000000000000000000000000ffffffffffffffffffffffff
    t1 = pow(c, 2, p)
    t1 = (t1 * c) % p
    t2 = pow(t1, 2**2, p)
    t2 = (t2 * t1) % p
    t3 = pow(t2, 2**4, p)
    t3 = (t3 * t2) % p
    t4 = pow(t3, 2**8, p)
    t4 = (t4 * t3) % p
    r = pow(t4, 2**16, p)
    r = (r * t4) % p
    r = pow(r, 2**32, p)
    r = (r * c) % p
    r = pow(r, 2**96, p)
    r = (r * c) % p
    return pow(r, 2**94, p)

def FindPointOnCurve(x):
    # This is the equation: y^2 = x^3-ax+b
    y2 = (x * x * x) - (3 * x) + P256.b
    y2 = y2 % P256.p
    y = p256_mod_sqrt(y2)
    return y2 == (y * y) % P256.p, y


def GeneratePrediction(observed, Q, d):
    checkbits = observed & 0xffff

    for high_bits in range(2**16):
        guess = (high_bits << (8 * 30)) | (observed >> (8 * 2))
        on_curve, y = FindPointOnCurve(guess)

        if on_curve:
            # use the backdoor to guess the next 30 bytes
            # point = Point(p256.curve, guess, y)
            point = Point(guess, y, curve=P256)
            s = (d * point).x
            r = (s * Q).x & (2**(8 * 30) - 1)

            # check the first 2 bytes against the observed bytes
            if checkbits == (r >> (8 * 28)):
                # if we have a match then we know the next 28 bits
                return r & (2**(8 * 28) - 1)

    return 0

main.py

seed = 0x1fc95c3714652fe2
P,Q,d = GenerateBackdoorConstants()
print("Generated these constants ", P, Q, d)
dualec = DualEC(seed, P, Q)

rngBits0 = dualec.GenerateBits()
rngBits1 = dualec.GenerateBits()
observed = (rngBits0 << (2 * 8)) | (rngBits1 >> (28 * 8))
print('Observed 32 bytes:\n{:x}'.format(observed))

print("\n\nRunning Backdoor\n")
predicted = GeneratePrediction(observed, Q, d)
print('Predicted 28 bytes:\n{:x}'.format(predicted))

actual = rngBits1 & (2**(8 * 28) - 1)
print('Actual 28 bytes:\n{:x}'.format(actual))