FullRepresentInteger.py

"""
An implementation of FullRepresentInteger from New algorithms for the Deuring correspondence.
https://eprint.iacr.org/2022/234

Provides gamma in O0 (instead of just ZZ[i,j]). This slows down Com and DecVrfy, but the corresponding isogegies are indistinguishable from random.

As an implementation note, we try and append "Heuristic" to function 
names when we only know solutions can be derived heuristically. If no 
solution is found, `None` is returned.
"""

# Sage imports
from sage.all import (
    ZZ,
    floor,
    sqrt,
    randint
)

# Local imports
from ideals import (
    quadratic_norm
)
from utilities import Cornacchia

from setup import (
    O0,
    p,
    q,
    j,
    ω
)

def RepresentIntegerHeuristic(M, parity=False):
    """
    Given an integer M, with M > p, attempts to
    find a random element γ with norm M.

    If no element is found after `bound` tries,
    returns none
    """

    def RepresentInteger(M, z, t, parity=False):
        M_prime = 4 * M - p * quadratic_norm(z, t)
        two_squares = Cornacchia(M_prime, -ZZ(ω**2))
        if two_squares:
            x, y = two_squares
            if parity and (x + t) % 2 == 1 and (y + z) % 2 == 1:
                return None
            return (x + ω * y + j * (z - ω * t)) / 2
        # No solution for the given M
        return None

    if M <= p:
        raise ValueError(f"Can only represent integers M > p.")
    m = max(floor(sqrt(4 * M / p)), 5)

    # TODO: how many times should we try?
    for _ in range(m**2):
        z = randint(-m, m)
        m2 = max(floor(sqrt(4 * M / p - z**2)), 5)
        t = randint(-m2, m2)
        γ = RepresentInteger(M, z, t, parity=parity)

        if γ is not None:
            # Found a valid solution, return
            assert γ.reduced_norm() == M, "The norm is incorrect"
            assert γ in O0, "The element is not contained in O0"
            return γ

    # No solution found, return None
    print(f"DEBUG [RepresentIntegerHeuristic]: No solution found")
    return None