We played Cyber Apocalypse 2021 and I have attempted several crypto challenges. I'll include the challenges Wii Phit and Hyper Metroid in this writeup.

Wii Phit

The aliens have encrypted our save file from Wii Phit and we’re about to lose our 4,869 day streak!! They’re even taunting us with a hint. I think the alien’s are getting a bit over-confident if you ask me.

CHALLENGE
from Crypto.Util.number import bytes_to_long
from secrets import FLAG,p,q

N = p**3 * q
e = 0x10001
c = pow(bytes_to_long(FLAG),e,N)

print(f'Flag: {hex(c)}')

# Hint

w = 25965460884749769384351428855708318685345170011800821829011862918688758545847199832834284337871947234627057905530743956554688825819477516944610078633662855
x = p + 1328
y = p + 1329
z = q - 1

assert w*(x*z + y*z - x*y) == 4*x*y*z

Challenge Analysis

Well... The code is pretty simple. It is just importing the number $p$ and $q$, and perform RSA encryption on the flag. Although the public key $N \neq pq$ as textbook RSA, it won’t be a big problem as we could still find it’s totient function $\varphi$ easily and calculate it’s private key $d$ if we could find $p$ and $q$.

So here’s the question, how could we find $p$ and $q$? Let’s see the hint.

Finding $p$ and $q$

We could see that the hint gives us four variables, $w$, $x$, $y$ and $z$, where all of them are either in terms of $p$, $q$ and some constant. The program then asserts that:

\[w(xz + yz - xy) = 4xyz\]

So lets expand it and make it become a multivariable equation in terms of $p$ and $q$ (and constant $w$):

\[ \begin{aligned} f(p, q) = &-4p^2q - (w-4)p^2 + (2w-10628)pq\\ &\quad-(2659w-10628)p + (2657w-7059648)q - (1767569w-7059648)\\ =&\ 0 \end{aligned} \]

Well it seems so complicated... To be honest, I had no idea what I should do next when I expanded the equation like this. I had tried various methods such as making it into modular equations. But of course as expected, it didn’t work.

When I was struggling with this equation, I suddenly noticed that there are only 2 positive terms in this equation, while the remaining 4 terms are negative terms. More importantly, by observing the degree of each terms, it seems that it's always likely that $f(p, q)$ will be negative no matter what $p$ and $q$ are. Therefore, I suspect that the maximum value of $f(p, q)$ would be 0 actually.

So lets assume $f(p_0, q_0)$ be the maximum point of $f(p, q)$. Then we can do partial differentiation on $f(p, q)$ to find out $p_0$ and $q_0$.

\[ \begin{gathered} \tfrac{\partial f}{\partial p} =-8pq - (2w-8)p + (2w-10628)q-(2659w-10628)\\ \tfrac{\partial f}{\partial q} =-4p^2 + (2w-10628)p + (2657w-7059648) \end{gathered} \]

For maximum point, $f_p = 0$ and $f_q = 0$. It's difficult to solve $f_p$ as it's still having 2 variables. However, $f_q$ now remains $p$ as the only term as a quadratic equation! So we are able to solve for $p_0$!

Does the quadratic equation has real positive integer root? Yes! Here's the value of $p_0$:

p_0 = 12982730442374884692175714427854159342672585005900410914505931459344379272923599916417142168935973617313528952765371978277344412909738758472305039316830099

So let's plug the value of $p_0$ to $f_p = 0$, then we can solve for $q_0$:

q_0 = 4376511920801673769046982367789644084746600661635151104602579081967083768976309788885633491753761209012042953502416064276555378570438196809829053232168930363213412874907199642703512833211163084612185095343783021125928854760406110492494016250237689683218940269389627326164130063600050024089561671913076715913062539181517022517918888557966172136200366952844819842691032260278270002343732370611629511112677638724984890716752906711541254394792207074699964511035503381

Finally, we use the value of $p_0$ and $q_0$ to check whether $w(xz + yz - xy) = 4xyz$ and it fits the equation! So the maximum point of $f(p, q)$ is really 0 and we have successfully retrive $p$ and $q$.

Then the rest is easy. For textbook RSA, the totient function $\varphi(pq) = (p-1)(q-1)$. In this case, the public key is not $pq$ but $p^3q$, so the totient function would be

\[\varphi(p^3q) = p^2(p-1)(q-1).\]

Then we can get the private key $d$ like textbook RSA. Using the private key $d$, we could decrypt the flag. Flag: CHTB{Erdos-Straus-Conjecture}

Me: Erdos-Straus-Conjecture? What is that?

Hyper Metroid

Dropping a morph ball bomb, Samus cracked open the floor and dropped down into the guts of Phaaze. At the end of the tunnel is a locked chest containing the hyper beam upgrade. Samus found the encrypted key preserved in a ball of glowing biomass, but can’t decode it. Help Samus capture the flag so she can eradicate the alien invasion once and for all.

CHALLENGE
from secrets import flag

def alien_prime(a):
    p = (a^5 - 1) // (a - 1)
    assert is_prime(p)
    return p


def encrypt_flag():
    e = 2873198723981729878912739
    Px = int.from_bytes(flag, 'big')
    P = C.lift_x(Px)
    JP = J(P)
    return e * JP


def transmit_point(P):
    mumford_x = P[0].list()
    mumford_y = P[1].list()
    return (mumford_x, mumford_y)


a = 1152921504606846997
alpha = 1532495540865888942099710761600010701873734514703868973
p = alien_prime(a)

FF = FiniteField(p)
R.<x> = PolynomialRing(FF)

h = 1
f = alpha*x^5

C = HyperellipticCurve(f,h,'u,v')
J = C.jacobian()
J = J(J.base_ring())

enc_flag = encrypt_flag()

print(f'Encrypted flag: {transmit_point(enc_flag)}')

#Encrypted flag: ([1276176453394706789434191960452761709509855370032312388696448886635083641, 989985690717445420998028698274140944147124715646744049560278470410306181, 1], [617662980003970124116899302233508481684830798429115930236899695789143420, 429111447857534151381555500502858912072308212835753316491912322925110307])

The missing piece

In the challenge, we are given a hyperelliptic curve $\mathcal{C}$, a constant $e$ and a point $Q$ (the encrypted flag). The objective is to find the point $P$ such that $Q = eP$. Similar to elliptic curves and modular exponentation, we are able to compute the inverse of $e$ to compute $P$ once we have the order of the field.

The missing piece of this puzzle is the order of the Jacobian, $\#J(\mathcal{C})$ and there is no efficient algorithm to calculate it generally. That means, this hyperelliptic curve $\mathcal{C}$ must be some special case so that we could caluclate its $\#J(\mathcal{C})$ easily.

Using my superior OSINT skills, I was able to find this special case in section 4.3.1 in A Comparison of Point Counting methods for Hyperelliptic Curves over Prime Fields and Fields of Characteristic 21.

As the same case mentioned in the paper, the prime $p$ generated in this challenge is a generalised Mersenne Prime (in a form of $\frac{a^5-1}{a-1}$) and $\alpha$ is also chosen such that

\[\alpha^{(p-1)/5} \equiv a\ (\text{mod}\ p).\]

This is exactly what we want and we could find out $\#J(\mathcal{C})$ in this challenge if we follow the paper!

But wait... multiplicative map?... automorphism?... What is $\zeta$? Oh, it's roots of unity, ...what is that?

I don't have much knowledge on mathematics that I totally can't understand what's the paper saying. Therefore, I found Mystiz in order to help but it seems that he can't understand the exact steps too... We then stucked there for about one and a half day just because we don't understand the mathematics in the paper.

Cryptosint

In the paper, the author gave an example of the Mersenne Prime method and calculated the $\#J(\mathcal{C})$ in that example. However, it doesn't show the steps in the example. By skimming the paper, I suddenly noticed that at the end of the paper, the author stated that all code used for the paper is available at the website2. Unfortunately, the paper is published in 2004 and the webpage has disappeared already.

Using my superior OSINT skills again. As I think that the website is an university website, it's having a high probability that it's being archived on the Internet. Therefore, I used Wayback Machine and I successfully found that webpage3 and I could even download the source code!

I then open the C source code on calculating the Mersenne Prime method and changed it to sage.

def calculateJacobiSum(a, genus, n):
    kth_root = 0
    index = 0
    temp_array = [None] * n
    jacobi = [None] * n
    
    """
    Evaluate the +/- \sigma^k part of the Jacobi sum
    This is only for genus 2 and 3 curves
    Maybe replace this by a formula later
    """
    modu = a % n
    if genus == 2:
        if modu == 0: kth_root = -1
        elif modu == 2: kth_root = -4
        elif modu == 3: kth_root = 2
        elif modu == 4: kth_root = 3
        else: kth_root = 0
    else:
        if modu == 0: kth_root = 4
        elif modu == 2: kth_root = -3
        elif modu == 3: kth_root = -5
        elif modu == 4: kth_root = 1
        elif modu == 5: kth_root = -6
        elif modu == 6: kth_root = 2
        else: kth_root = 0
    
    #Calculate rest of jacobi sum
    if genus == 2:
        temp_array[0] = a*a
        temp_array[1] = -a
        temp_array[2] = 0
        temp_array[3] = -a
        temp_array[4] = 1
    else:
        temp_array[0] = a*a*a
        temp_array[1] = -a*a
        temp_array[2] = a
        temp_array[3] = -1
        temp_array[4] = -a*a
        temp_array[5] = a - (a*a)
        temp_array[6] = a
    
    #Now multiply by kth_root
    for i in range(n):
        index = (abs(kth_root) + i) % n
        jacobi[index] = temp_array[i]
        if kth_root < 0:
            jacobi[index] *= -1
    
    return jacobi

def evaluateOrderTwists(jacobi, i, j, n):
    twist = 0
    
    order = [None] * (n - 1)
    roots_unity = [None] * n
    
    #Calculate roots of unity that are needed
    roots_unity[0] = 1
    for k in range(1, n):
        real_part = cos((2*k*pi)/n)
        imag_part = sin((2*k*pi)/n)
        roots_unity[k] = real_part + imag_part*I
    
    #Add in the "twist" to the jacobi sum
    twist = pow(-1, i)
    jacobi[j] = jacobi[j] + twist;
    
    print("twist: " + str(twist))
    
    order[0] = 0
    #Evaluate first iteration
    for k in range(n):
        order[0] += (jacobi[k] * roots_unity[k])
    
    #Evaluate norm of jacobi sum (and twist)
    for k in range(2, n):
        jac2 = [None] * n
        index = 0
        jac2[0] = jacobi[0]
        
        #Bump all the roots of unity up by the "k-th" power
        for l in range(1, n):
            index = (l * k) % n
            jac2[index] = jacobi[l]

        order[k-1] = 0
        #Evaluate k-th iteration
        for l in range(n):
            order[k-1] += (jac2[l] * roots_unity[l])
    
    #Find the total order of the jacobian
    total_order = order[0] * order[1] * order[2] * order[3]

    return total_order.real()

You ask me that is the i and j in evaluateOrderTwists? I don't know either. Maybe let's just follow the input of the example in the webpage.

Therefore, by running the code above in sage, I've successfully found the order of the jacobian.

#J(C) = 3121748550315992691688089417819100427615767189940196696753029154066816327228360361217733781238925753559457441594415001909435541252098381221684880

Then the remaining thing is lift the encrypted flag point onto the jacobian, and multiply it by the inverse to get back the original flag point. Then we could get the flag.

Here's the solve script:

def alien_prime(a):
	p = (a^5 - 1) // (a - 1)
	assert is_prime(p)
	return p

def calculateJacobiSum(a, genus, n):
    kth_root = 0
    index = 0
    temp_array = [None] * n
    jacobi = [None] * n
    
    """
    Evaluate the +/- \sigma^k part of the Jacobi sum
    This is only for genus 2 and 3 curves
    Maybe replace this by a formula later
    """
    modu = a % n
    if genus == 2:
        if modu == 0: kth_root = -1
        elif modu == 2: kth_root = -4
        elif modu == 3: kth_root = 2
        elif modu == 4: kth_root = 3
        else: kth_root = 0
    else:
        if modu == 0: kth_root = 4
        elif modu == 2: kth_root = -3
        elif modu == 3: kth_root = -5
        elif modu == 4: kth_root = 1
        elif modu == 5: kth_root = -6
        elif modu == 6: kth_root = 2
        else: kth_root = 0
    
    #Calculate rest of jacobi sum
    if genus == 2:
        temp_array[0] = a*a
        temp_array[1] = -a
        temp_array[2] = 0
        temp_array[3] = -a
        temp_array[4] = 1
    else:
        temp_array[0] = a*a*a
        temp_array[1] = -a*a
        temp_array[2] = a
        temp_array[3] = -1
        temp_array[4] = -a*a
        temp_array[5] = a - (a*a)
        temp_array[6] = a
    
    #Now multiply by kth_root
    for i in range(n):
        index = (abs(kth_root) + i) % n
        jacobi[index] = temp_array[i]
        if kth_root < 0:
            jacobi[index] *= -1
    
    return jacobi

def evaluateOrderTwists(jacobi, i, j, n):
    twist = 0
    
    order = [None] * (n - 1)
    roots_unity = [None] * n
    
    #Calculate roots of unity that are needed
    roots_unity[0] = 1
    for k in range(1, n):
        real_part = cos((2*k*pi)/n)
        imag_part = sin((2*k*pi)/n)
        roots_unity[k] = real_part + imag_part*I
    
    #Add in the "twist" to the jacobi sum
    twist = pow(-1, i)
    jacobi[j] = jacobi[j] + twist;
    
    print("twist: " + str(twist))
    
    order[0] = 0
    #Evaluate first iteration
    for k in range(n):
        order[0] += (jacobi[k] * roots_unity[k])
    
    #Evaluate norm of jacobi sum (and twist)
    for k in range(2, n):
        jac2 = [None] * n
        index = 0
        jac2[0] = jacobi[0]
        
        #Bump all the roots of unity up by the "k-th" power
        for l in range(1, n):
            index = (l * k) % n
            jac2[index] = jacobi[l]

        order[k-1] = 0
        #Evaluate k-th iteration
        for l in range(n):
            order[k-1] += (jac2[l] * roots_unity[l])
    
    #Find the total order of the jacobian
    total_order = order[0] * order[1] * order[2] * order[3]

    return total_order.real()

def data_to_jacobian(data):
    u, v = data[0], data[1]
    u = sum([u[i] * x^i for i in range(3)])
    v = sum([v[i] * x^i for i in range(2)])

    return J((u, v))

a = 1152921504606846997
g = 2
n = 5
e = 2873198723981729878912739

jacobi = calculateJacobiSum(a, g, n)

#order = evaluateOrderTwists(jacobi, 0, 1, n).n(digits=200)
order = 3121748550315992691688089417819100427615767189940196696753029154066816327228360361217733781238925753559457441594415001909435541252098381221684880

alpha = 1532495540865888942099710761600010701873734514703868973
p = alien_prime(a)

FF = FiniteField(p)
R.<x> = PolynomialRing(FF)

C = HyperellipticCurve(alpha*x^5,1,'u,v')
J = C.jacobian()
J = J(J.base_ring())

enc_flag = ([1276176453394706789434191960452761709509855370032312388696448886635083641, 989985690717445420998028698274140944147124715646744049560278470410306181, 1], [617662980003970124116899302233508481684830798429115930236899695789143420, 429111447857534151381555500502858912072308212835753316491912322925110307])

inv = inverse_mod(e, order)
J_point = data_to_jacobian(enc_flag)
flag_point = inv*J_point
flag_int = flag_point[0].coefficients()[0]

print(-flag_int % p)

Flag: CHTB{hyp3r_sp33d_c0unting!!}

Who need maths if you can OSINT?


  1. Colm Ó hÉigeartaigh (2004) "A Comparison of Point Counting methods for Hyperelliptic Curves over Prime Fields and Fields of Characteristic 2"
    https://eprint.iacr.org/2004/241.pdf [return]
  2. Colm Ó hÉigeartaigh (2004) "Hyperelliptic Curve Cryptographic Software" (Invalid)
    http://www.computing.dcu.ie/~coheigeartaigh/crypto.html [return]
  3. Colm Ó hÉigeartaigh (2004) "Hyperelliptic Curve Cryptographic Software"
    https://web.archive.org/web/20050924125637/http://www.computing.dcu.ie/~coheigeartaigh/crypto.html [return]