# BarakLe challenge revient à résoudre une instance de ECDLP sur une courbe Hessienne.
La courbe est d'équation$$    (E): x^3 + y^3 + c - dxy = 0$$avec```pythonp = 73997272456239171124655017039956026551127725934222347d = 68212800478915688445169020404812347140341674954375635c = 1```
Note, ce n'est techniquement **pas** une courbe elliptique à cause du $y^3$.
Pour les curieux, voici à quoi ressemble une courbe Hessienne:![E](https://i.imgur.com/WMg1Ozo.png)C'est très similaire au Folium de Descartes, mais c'est pas exactement la même courbe.
![Point addition](https://i.imgur.com/wOnowOh.png)L'addition de points est similaire à l'addition classique pour les courbes elliptiques.On trace la droite passant par les deux points, on trouve la 3e intersection eton effectue une symétrie axiale, cette fois-ci autour de $(\Delta): y=x$.
Il nous est donné deux points $P$, $Q$ de sorte que $P = Q\cdot m$ pour un certain $m$.L'objectif est de retrouver $m$, qui encode la solution du challenge.
```P = (71451574057642615329217496196104648829170714086074852,     69505051165402823276701818777271117086632959198597714)Q = (40867727924496334272422180051448163594354522440089644,     56052452825146620306694006054673427761687498088402245)```
## Trouver l'ordre du point P
Il est possible de construire une équivalence birationnelle entre $(E)$ et unecourbe elliptique[1].Le script suivant illustre la construction de cette courbe elliptique:
```python=p = 73997272456239171124655017039956026551127725934222347G = GF(p)d = G(68212800478915688445169020404812347140341674954375635)D = d/3c = 1# Hessian Curve:# x³ + y³ + 1 = 3Dxy
# Weierstrass equivalenta = -27*D*(D^3 + 8)b = 54*(D^6 - 20*D^3 - 8)E2 = EllipticCurve(G, [a, b])print(E2)# Elliptic Curve defined by y^2 = x^3 + 44905983883632632311912975168565494049729462391119290*x + 4053170785171018449128623853386306889464200866918538 over Finite Field of size 73997272456239171124655017039956026551127725934222347```
Les deux courbes étant équivalentes, cela nous assure qu'elles possèdent le mêmeordre, que nous pouvons calculer grâce à sage:
```python=order = E2.order()# 73997272456239171124655016995459084401465136460086688```
Par théorème de Lagrange, l'ordre de $P$ doit obligatoirement diviser l'ordre de $E$,donc l'ordre de $E2$.Son ordre va donc être le plus petit diviseur $k$ de `order` tel que $P\cdot k = O$.
Le script suivant implémente les divers opérateurs sur la courbe Hessienne,puis calcule l'ordre $n$ du point $P$:
```python=class HessianPoint:    def __init__(self, x, y):        self.x = G(x)        self.y = G(y)
        assert (x, y) == (0, 0) or x^3 + y^3 + c - d*x*y == 0
    def __eq__(P, Q):        return (P.x, P.y) == (Q.x, Q.y)
    def __add__(P, Q):        if (P.x, P.y) == (0, 0):            return Q        if (Q.x, Q.y) == (0, 0):            return P        if (P.x, P.y) == (Q.y, Q.x):            return HessianPoint(0, 0)
        x1, y1 = P.x, P.y        x2, y2 = Q.x, Q.y        if P == Q:            x3 = y1 * (c - x1^3) / (x1^3 - y1^3)            y3 = x1 * (y1^3 - c) / (x1^3 - y1^3)        else:            x3 = (y1^2*x2 - y2^2*x1) / (x2*y2 - x1*y1)            y3 = (x1^2*y2 - x2^2*y1) / (x2*y2 - x1*y1)        return HessianPoint(x3, y3)
    def __mul__(P, m):        if (P.x, P.y) == (0, 0):            return P        R = HessianPoint(0, 0)        while m != 0:            if m & 1:                R = P + R            m = m >> 1            if m != 0:                P = P + P        return R
    def __neg__(self):        return HessianPoint(self.y, self.x)
    def __rmul__(P, m):        return P*m
    def __repr__(self):        return f"({self.x}, {self.y})"
    def __hash__(self):        return hash((self.x, self.y)) P = HessianPoint(71451574057642615329217496196104648829170714086074852, 69505051165402823276701818777271117086632959198597714)Q = HessianPoint(40867727924496334272422180051448163594354522440089644, 56052452825146620306694006054673427761687498088402245)O = HessianPoint(0, 0)
for k in order.divisors():    if P * k == O:        breakn = k# 3083219685676632130193959041477461850061047352503612```
On remarque que l'ordre de $P$ se factorise simplement:$$    n = 2^2 × 3^2 × 17 × 2341 × 23497 × 500369 × 5 867327 × 33 510311 × 13824 276503 × 67342 255597$$
Le plus grand facteur premier de $n$ étant relativement petit, on peut se permettre de résoudre le logarithme discret avec une variante simplifiée de l'algorithme de Pohlig-Hellman.
## Résoudre le Logarithme Discret
Pohlig-Hellman nous permet de simplifier ECDLP en plusieurs plus petites instances pour chaque facteur de $n$ que l'on peut par la suite combiner grâce à un CRT.Voici une variante simplifiée de l'algorithme de Pohlig-Hellman qui nous est suffisant au vu des contraintes:
Pour chaque facteur $p_i^{e_i}$ de $n$, on cherche à trouver $k_i \equiv m \mod p_i^{e_i}$.Pour ce faire, on considère les deux points $P' = P \cdot \dfrac n {p_i^{e_i}}$ et $Q' = Q \cdot \dfrac n {p_i^{e_i}}$.Étant donné que $P' \cdot p_i^{e_i} = P \cdot n$, il est clair que l'ordre de $P'$ est de $p_i^{e_i}$.
Ainsi, on peut résoudre $P' \cdot k_i = Q'$ et obtenir la valeur de $m \mod p_i^{e_i}$.
Afin de résoudre cette plus petite instance, j'utilise l'algorithme baby-step giant-step,de complexité temporelle $O(\sqrt{p_i^{e_i}})$.
```python=def baby_step_giant_step(P, Q, order):    """    Résoud P·x = Q pour x entre 0 et order-1    """
    m = isqrt(order)    #   x = am + b    # ⇒ Q = P*(am + b)    # ⇒ P*(am) = Q+P*(-b)
    memory = {}    # On y stocke toutes les valeurs possibles de P*(am),    # associées à leur valeur de `a` respective.    M = O    S = P*m # giant step
    for a in range(m):        # M = P*(am)        memory[M] = a        M += S
    b = 0    M = Q    S = -P # baby step
    while M not in memory:        # M = Q+P*(-b)        b += 1        M += S
    a = memory[M]    x = a*m + b    assert P * x == Q    return x```
Et ici se trouve l'algorithme complet en action:
```python=factors = n.factor()p = [factor[0] for factor in factors]e = [factor[1] for factor in factors]
crt_a = []crt_b = []for i in range(len(p)):    pi = p[i]    ei = e[i]    P0 = (n // (pi^ei)) * P    Q0 = (n // (pi^ei)) * Q    k = baby_step_giant_step(P0, Q0, pi^ei)    crt_a.append(k)    crt_b.append(pi^ei)    print(f"m = {k} mod {pi^ei}")# m = 1 mod 4# m = 1 mod 9# m = 13 mod 17# m = 1489 mod 2341# m = 18879 mod 23497# m = 339241 mod 500369# m = 4831116 mod 5867327# m = 31953114 mod 33510311# m = 10804270204 mod 13824276503# m = 57315765050 mod 67342255597
m = CRT(crt_a, crt_b)# 1780694557271320552511299360138314441283923223949197assert P * m == Q```
## Retrouver le flag
Dernière subtilité, il existe plusieurs valeurs de $m$ valides tel que $P \cdot m = Q$.En effet $P \cdot (m + \lambda n) = P \cdot m + P \cdot n \cdot \lambda = Q + O \cdot \lambda = Q$ pour toute valeur de $\lambda$ entière.La seule restriction donnée dans le fichier d'entrée est que $m < p$, alors on testetoutes les valeurs de $m$ valides jusqu'à trouver le flag:
```python=p = 73997272456239171124655017039956026551127725934222347while m < p:    flag = long_to_bytes(m)    try:        print('CCTF{' + flag.decode() + '}')    except ValueError:        pass    m += n
# CCTF{_hE5S!4n_f0rM_0F_3CC!!}```
## Références
[1] https://link.springer.com/content/pdf/10.1007/3-540-44709-1_33.pdf
