Sage Code for GGH Cryptanalysis by Hu and Jia

Recently, Yupu Hu and Huiwen Jia put a paper on the Cryptology ePrint Archive which describes a successful attack of the GGH (and GGHLite) candidate multilinear map. The attack does not try to recover the secret g or any other secret parameter of the map. Instead, it solves the Extraction \kappa-graded CDH (Ext-GCDH) problem directly.

For GGHLite the Ext-GCDH problem is defined as computing the higher order bits of v = h/g \cdot \prod_{i=0}^{\kappa} e_i from u_i = e_i \cdot y + \rho_{i0} x_0 + \rho_{i1} x_1, x_0 = b_0 \cdot g/z, x_1 = b_1 \cdot g/z, y = (1+a \cdot g)/z and p_{zt} = h \cdot z^\kappa/g.

All elements involving a z live in R_q = \mathbb{Z}_q[x]/(x^n+1) where n is a power of two, the other elements live in R = \mathbb{Z}[x]/(x^n+1). All elements except for z and h are small: z is random in R_q and \Vert h \Vert \approx \sqrt{q}. We define “levels” by the powers of z, e.g. (e_0 +cg)/z^i is a level-i encoding of e_0 for some c \in R.

In a non-interactive key exchange (NIKE) based on GGH-like maps, the higher order bits of the value v are essentially the shared secret. Let’s look at a simple Sage implementation of non-interactive key echange based on a GGH-like multilinear map.

class NIKESloppy:
    """
    A sloppy implementation of NIKE for testing.

    By "sloppy" we mean that we do not care about distributions,
    verifying some properties or exact sizes.
    """

    def __init__(self, n, kappa):
        """
        Initialise a new sloppy GGHLite instance.

        :param n: dimension (must be power of two)
        :param kappa: multilinearity level `κ`
        """
        self.n = ZZ(n)
        if not n.is_power_of(2):
            raise ValueError
        self.kappa = ZZ(kappa)

        self.R = ZZ["x"]
        self.K = CyclotomicField(2*n, "x")
        self.f = self.R.gen()**self.n + 1

        self.sigma = self.n
        self.sigma_p = round(self.n**(ZZ(5)/2) * self.sigma)
        self.sigma_s = round(self.n**(ZZ(3)/2) * self.sigma_p**2)

        self.q = next_prime((3*self.n**(ZZ(3)/2) * self.sigma_s*self.sigma_p)**(3*self.kappa))
        self.Rq = self.R.change_ring(IntegerModRing(self.q))

        g = self.discrete_gaussian(self.sigma)
        # we enforce primality to ensure we can invert mod ⟨g⟩
        while not self.K(g).is_prime():
            g = self.discrete_gaussian(self.sigma)

        h = self.discrete_gaussian(sqrt(self.q))
        z = self.Rq.random_element(degree=n-1)

        # encoding of one
        a = self.discrete_gaussian(self.sigma_p/self.sigma)
        self.y = ((1+a*g) * z.inverse_mod(self.f)) % self.f

        # encoding of zero
        b0 = self.discrete_gaussian(self.sigma_p/self.sigma)
        self.x0 = ((b0*g) * z.inverse_mod(self.f)) % self.f

        # encoding or zero
        b1 = self.discrete_gaussian(self.sigma_p/self.sigma)
        self.x1 = ((b1*g) * z.inverse_mod(self.f)) % self.f

        # zero testing parameter
        self.pzt = h*z**kappa * self.Rq(g).inverse_mod(self.f)

        # we migh as well give these to the attacker, as they can be computed
        # easily
        self.G = self.canonical_basis(g)
        self.I = self.K.ideal([self.R(row.list()) for row in self.G.rows()])

Logarithms of Euclidean norms are useful as “small” and “not-so-small” are what distinguishes between secret and not-so-secret.

    @staticmethod
    def norm(f):
        "Return log of Euclidean norm of `f`"
        return log(max(map(abs, f)), 2.).n()

It is easy to compute a canonical basis for \langle g \rangle \in R from publically available information. In this simplified version we simply give this canonical basis to the attacker.

    def canonical_basis(self, g):
        "Return HNF of g"
        x = self.R.gen()
        n = self.n
        G = [x**n + ((x**i * g) % self.f) for i in range(n)]
        G = [r.list()[:-1] for r in G]
        G = matrix(ZZ, n, n, G)
        return G.hermite_form()

It will be useful below to compute short remainders modulo some number ring element. We use Babai’s rounding trick for that.

    def rem(self, h, g):
        """
        Compute a small remainder of `h` modulo `g`

        :param h: a polynomial in ZZ[x]/x^n+1
        :param g: a polynomial in ZZ[x]/x^n+1
        :returns: small representative of `h mod g`

        """
        g_inv = g.change_ring(QQ).inverse_mod(self.f)
        q = (h * g_inv) % self.f
        q = self.R([round(c) for c in q.coefficients()])
        r = (h - q*g) % self.f
        return r

We do not care about distributions of elements in our toy example, so we simply return random polynomials where the coefficients are bounded by \sigma.

    def discrete_gaussian(self, sigma):
        """
        sample element with coefficients bounded by sigma

        .. note:: this is obviously *not* a discrete Gaussian, but we don't care here.
        """
        sigma = round(sigma)
        return self.R.random_element(x=-sigma, y=sigma, degree=self.n-1)

We need a function to sample new encodings:

    def sample(self):
        """
        Sample level-0 and level-1 encodeing.
        """
        e0 = self.discrete_gaussian(self.sigma_p)
        r0 = self.discrete_gaussian(self.sigma_s)
        r1 = self.discrete_gaussian(self.sigma_s)

        u0 = ((e0*self.y + r0*self.x0 + r1*self.x1)) % self.f
        return e0, u0

For the final step we need to extract a representation over the Integers from our representation modulo q. Sage does not normalise to [-\lfloor q/2 \rfloor,\dots, \lfloor q/2 \rfloor] but to [0,\dots,q). We work around this behaviour. Also, lower order bits are noisy so we discard them.

    def extract(self, p, mul_pzt=True, drop_noise=True):
        """
        Convert mod q element to element over the integers

        :param p: a polynomial in ZZ[x]/x^n+1
        :param mul_pzt: multiply by pzt
        :param drop_noise: remove lower order terms
        """
        p = p % self.f
        if mul_pzt:
            p = (self.pzt*p) % self.f
        q = parent(p).base_ring().order()
        f = []
        for c in p.change_ring(ZZ).list():
            if c<q/2:
                f.append(c)
            else:
                f.append(c-q)

        if drop_noise:
            f = [f[i] >> floor(log(self.q//2, 2)) for i in range(self.n)]
        return ZZ['x'](f)

Finally, we implement NIKE: we multiply \kappa level-1 encodings and our secret level-0 encoding to produce a level-\kappa encoding of \prod_{i=0}^\kappa e_i from which we then extract our shared key.

    def __call__(self, e0, U):
        """
        Run NIKE with `e_0` and `U`

        :param e0: a level-0 encoding
        :param U: κ level-1 encodings
        """
        if len(U) != self.kappa:
            raise ValueError

        t = prod(U) % self.f
        t = t * e0 % self.f
        t = self.extract(t)
        return tuple(t.list())

Let me stress once more that the above is not a proper implementation of a GGH-like candidate multilinear maps. For starters, the method called discrete_gaussian does not return elements following a discrete Gaussian distribution. We also skipped several tests from the GGH and GGHLite specifications. However, this simplified variant is sufficient to explain the attack of Hu and Jia.

The attack proceeds in two steps. The first step was already known as the “weak discrete log attack”. (Update: Hu and Jia do not present the first step of their attack this way. Instead, they give another way of computing a “kind of” level-0 encoding from a, say, level-1 encoding. Damien Stéhle was the first to point out to me that this step can be replaced by the “normal” weak discrete log attack.) Say, we have some level-1 encoding u_0 = e_0 \cdot y + \rho_{00} x_0 + \rho_{01} x_1. We can compute:

v_0 = \left(u_0 \cdot x_0 \cdot y^{\kappa-2} \cdot p_{zt}\right) \mod q = e_0\, b_0\, h + \xi_0 g

for some \xi_0. Note that the right hand side has no modular reduction modulo q. Also note that e_0\, b_0\, h + \xi_0 g \equiv  e_0\, b_0\, h \mod \langle g \rangle, i.e. we have a representative of e_0\, b_0\, h.

Similarly, we can compute:

v_{-1} = \left(x_0 \cdot y^{\kappa-1} \cdot p_{zt}\right) \mod q = b_0\, h + \xi_1 g  \equiv  b_0\, h \mod \langle g \rangle

Now, assume v_{-1} is invertible modulo \langle g\rangle, then we can compute \tilde{v}_{-1} = v_{-1}^{-1} \mod \langle g \rangle and \tilde{v}_{-1} \cdot v_0 \mod \langle g \rangle \equiv e_0 \mod \langle g \rangle. In other words, we can compute a representative of e_0 “kind of” at level-0. Here’s some Sage code implementing this step:

def weak_dlog(params, u):
    """
    Weak discrete log attack from [EC:GarGenHal13]_

    :param params: GGH parameters
    :param u: a level-1 encoding
    """
    kappa = params.kappa
    v0 = params.extract(u*params.x0*params.y**(kappa-2), drop_noise=False)
    vn = params.extract(params.y**(kappa-1)*params.x0, drop_noise=False)

    v0 = params.K(v0).mod(params.I)
    vn = params.K(vn).mod(params.I)

    r =  (v0 * vn.inverse_mod(params.I)).mod(params.I)
    return r.polynomial().change_ring(ZZ)

Now, compute representivatives v_i for all i \in \{0,\dots,\kappa\} and finally \eta = \prod_{i=0}^{\kappa} v_i. This produces a “kind of” level-0 encoding of \prod_{i=0}^{\kappa} e_i, where “kind of” stands for “not small”. If \eta was small, we could simply compute \eta \cdot y^{\kappa} \cdot p_{zt} to solve the Ext-GCDH problem. However, \eta is not small.

This is where Hu and Jia come in with a clever and remarkably simple step two. They define:

  • Y   = \left(y^{\kappa-1} x_1     p_{zt}\right) \mod q = h \cdot (1+ag)^{\kappa-1}\cdot b_1
  • X_i = \left(y^{\kappa-2} x_i x_1 p_{zt}\right) \mod q = h \cdot (1+ag)^{\kappa-2}\cdot b_1 \cdot b_i \cdot g

Again, note that the right hand sides have no modular reductions modulo q because they are “somewhat short”.

Recall that \eta = \prod_{i=0}^{\kappa} e_i + \zeta_0 g for some \zeta_0 \in R and compute:

  1. \eta' := Y \cdot \eta= Y\cdot\prod_{i=0}^\kappa e_i + \zeta_1\, b_1\, g. Note that Y is a multiple of b_1\, g.
  2. \eta'' := \eta' \mod X_1 = Y\cdot\prod_{i=0}^\kappa e_i + \zeta_2\, b_1 g which holds because X_1 is a multiple of b_1\,g. Note that \eta'' has roughly the same size as X_1, which is “somewhat short”.
  3. \eta''' := \eta' \cdot y \cdot (x_1)^{-1} \eta'' \mod q This effectively replaces the factor x_1 in Y by y. Consequently, we have

    \eta''' = y^\kappa\cdot p_{zt} \cdot \prod_{i=0}^\kappa e_i + \zeta_2 (1+ag) \mod q

    where the second summand is small. In other words, we have computed

    h/g\cdot \prod_{i=0}^\kappa e_i + \textnormal{something small} \mod q,

    i.e. we have solved the Ext-GCDH problem.

Here’s some Sage code implementing these three steps:

def ggh_break(params, U):
    """
    Attack from [EPRINT:HuJia15]_

    :param params: GGH parameters
    :param U: κ+1 level-1 encodings
    """
    if len(U) != params.kappa + 1:
        raise ValueError

    kappa = params.kappa

    # large "level-0" encodings of U
    V = [weak_dlog(params, u) for u in U]

    # large "level-0" encoding of prod(U)
    eta = prod(V) % params.f

    # Y = (b_0 ⋅ h ⋅ (1+ag))^(κ-2) ⋅ (1+ag)
    Y = params.extract(params.y**(kappa-1) * params.x0, drop_noise=False)
    # X = (b_0 ⋅ h ⋅ (1+ag))^(κ-2) ⋅ (b_0⋅g)
    X = params.extract(params.y**(kappa-2) * params.x0**2, drop_noise=False)

    # η'' = (Y ∏ u_i) mod X
    eta_ = (eta * Y) % params.f
    eta__ = params.rem(eta_, X)

    # η''' = η'' ⋅ y/x
    xy = params.y * params.x0.inverse_mod(params.f) % params.f
    eta___ = params.extract(xy * eta__, mul_pzt=False)
    return tuple(eta___.list())

Finally, let’s put everything together and run the attack. If everything is well, then the following code should print the same tuple three times.

def run(n, kappa):
    params = NIKESloppy(n, kappa)
    EU = [params.sample() for i in range(kappa+1)]
    E = [eu[0] for eu in EU]
    U = [eu[1] for eu in EU]
    # fist make sure NIKE is correct, by running it twice
    print params(E[0], U[1:])
    print params(E[-1], U[:-1])
    # now run the attack
    print ggh_break(params, U)

One thought on “Sage Code for GGH Cryptanalysis by Hu and Jia

Leave a comment