Another year, another hxp CTF. Here's a writeup for the evidently hardest crypto challenge "f_cktoring", organically sourced directly from the author!

This challenge was a classical RSA decryption challenge: Given a public key and ciphertext, recover the message.

Because hxp 💝 you and guessing 🟰 stupid, we obviously provided the full source code of the SageMath script that generated the data..:

#!/usr/bin/env sage
R = __import__('random').SystemRandom().randint

def gen_prime():
    while True:
        p = (ZZ^2).0
        while p.norm() < 5^55:
            a,b = ((-1)^R(0,1)*R(1,7^y) for y in (2,1))
            p *= matrix([[a,b],[123*b,-a]])
        p += (ZZ^2).0
        p *= diagonal_matrix((1,123)) * p
        if is_pseudoprime(p): return p

n = prod(gen_prime() for _ in 'yyyyyyy')
print(f'{n = :#x}')

m = int.from_bytes(open('flag.txt','rb').read().strip(), 'big')
c = int(pow(m, 1|2^2^2^2, n))
print(f'{c = :#x}')

Output:

n = 0xb96069105fd32b73363cdf8eae91b1f5ec165b1ffaa8b9b56af7414e8cbeac03b897e8f8bea3067d7cd2c4baa4d1cedff94c6918fb266039057843f1ee43bd6e3b0fcb8eed1e0945a4eacb4f4c720a6fffaaadb14b6ef83549868db5a0fb39f29b97a94e1449f4b56a7775fc10cbe9e2a0d39c74ffcde60bf66f8b880654d4eb526283ab37a72e03fe1a68af0d2017ae04f7d3be6c2fad768d76a67ef2ceacac2b846ec3393f6c962861d06105134c5b7884bf561401a9fd897fff7aeb4c591a3cbd9f72b68eb4ae0ac64f213ef38449ca7f6ae98113a95731525661d5bab460b1b9f3eb801
c = 0x2c51316f8ace726cb13f27470717f85604187480bf3fbb7d5a7a2b117d7e93d982370a204705d66542d75dee2fe4e6ad1fbab2406dd7c8ca188553404f4ae7e713f18ae4aa1987d85d43e56e1ed15b9baddc24675c8db1fe876f7f579c73087ecf3d62359fabcfba7eb86e6e820f49791f5766edc6b519e1e402ccc252c2c70de3ad199b11227470fba03d34b5eb91a8b2af2f5f599e40cd1075ed5087395b6041fb720f08c43d508e493a48bb886f969d45f33324ad09d8bcb99abcef47339f6f34a38ed133277e3199935063c49141f9a5265ef3345ee95dbd896c473a09114913b968d54

What's going on?

Let's unpack what this code does:

• The final RSA modulus $n$ is a product of seven primes generated by gen_prime(), and everything else is just normal RSA.
• The main loop of gen_prime() starts off with a tuple $(x_0,y_0)=(1,0)$. It then repeatedly samples two integers $(a,b)$ of absolute values in the range $[1;49]$ resp. $[1;7]$ and replaces the tuple $(x_i,y_i)$ by the tuple \[ (x_{i+1},y_{i+1}):=(ax_i+123by_i, bx_i-ay_i) \,\text. \] This is repeated until the integers $(x,y)$ reach a length of about 128 bits. Finally, the function returns \[ p := (x{+}1)^2+123y^2 \] assuming this number is prime; else it simply tries again.

That's all.

Quadratic field trip!

One key ingredient to solving the challenge is to understand what the weird inner loop actually means. A hint towards this can be found in the shape of the final expression $(x{+}1)^2+123y^2$, which just happens to equal the norm of the algebraic integer $(x{+}1) + \sqrt{-123}\cdot y$!

Knowing this, it should not be unreasonably hard to recognize that the weird matrix multiplication in the inner loop simply computes a multiplication in the ring $\mathbb Z[\sqrt{-123}]$: Writing $\tau:=\sqrt{-123}$, we have \[ (x_i-\tau y_i)\cdot (a+\tau b) = x_ia + \tau (x_ib-y_ia) - \tau^2y_ib = (x_ia + 123y_ib) + \tau (x_ib-y_ia) = x_{i+1} + \tau \cdot y_{i+1} \,\text. \]

Now, since the norm is multiplicative and invariant under complex conjugation $\tau\mapsto-\tau$, one immediate consequence of this is that the norm $x^2+123y^2$ of the last tuple $(x,y)$ is the same as the product of all the norms of the elements $(a+\tau b)$ used in the loop. Among all the facts this implies, the most important one is that

(This means that all prime factors are "small", in this case upper bounded by $271723\approx 2^{18}$ and typically even smaller. This is because $a,b$ are pretty small, hence the norm of $a+\tau b$ is as well, and multiplicativity does the rest.)

Factoring using algebraic groups

To summarize, we now know that something with a very clear algebraic relation to our desired prime $p$ is a smooth integer: There exist $x,y\in\mathbb Z$ such that $p=(x{+}1)^2+123y^2$ whereas the integer $x^2+123y^2$ is known to be smooth.

Let's recap some well-known factorization algorithms which exploit smoothness. All of the following algorithms have one common underlying recipe:

1. Find an algebraic group (computed modulo $n$) whose reduction modulo special prime divisors $p\mid n$ has somewhat predictable group order.
2. Compute an exponentiation in that group to a power that will result in a special element modulo all the special prime divisors of $n$.
3. Use a $\gcd$ with $n$ involving the special element to split off the special factors from $n$.

The way the exponent for step 2 is chosen is as a likely multiple of (some particular divisor of) the (unknown) order of the base modulo a special prime factor. In practice, we hope that the order is smooth and let the exponent be a large smooth integer.

Note (let that sink in): The reason any of this has to do with smoothness is that the $\mathrm{lcm}$ of smooth integers is dramatically smaller than for random integers of comparable size. In other words, multiples of unknown smooth integers are relatively easy to guess.

"Basic" instantiations of this approach:

• Pollard's $p-1$ algorithm finds prime factors $p$ for which $p-1$ is smooth by computing exponentiations in the ring $\mathbb Z/n$. The hope is to successfully guess some multiple $e$ of $p-1$, which implies that the power $a^e$ will be congruent to $1$ modulo the (unknown) prime $p$, and thus we can split off $p$ from $n$ by computing $\gcd(a^e-1,n)$.
• Williams' $p+1$ algorithm finds prime factors $p$ for which $p+1$ is smooth by computing exponentiations in a degree-$2$ extension of $\mathbb Z/n$, say $R=\mathbb Z[T]/(n, T^2{-}A)$ for some integer $A$. There's a good chance (about half of all $A$) that the reduction of $R$ modulo $p$ is $\mathbb F_{p^2}$, whose multiplicative group has order $p^2-1=(p-1)(p+1)$. Suppose that we successfully guess a multiple $e$ of $p+1$: Picking some base element $a=a_0+a_1T$ at random and computing $b_0+b_1T:=a^e$ in $R$ then yields an element whose reduction modulo $p$ lies in the subgroup $\mathbb F_p^\times$ of $\mathbb F_{p^2}^\times$. In more concrete terms, this means that $b_1 \equiv 0 \pmod p$, and therefore computing $\gcd(b_1,n)$ splits off $p$.
• There exists a generalization of this method for higher degrees ($p-1$ was degree one, $p+1$ was degree two) due to Bach and Shallit. The technical condition turns into the smoothness of certain divisors of $p^k-1$, but the underlying idea is exactly the same as for $p-1$ and $p+1$, just with higher-degree polynomials.

All of the above have one big problem (from the factorer's perspective): You'd have to be astronomically lucky to be facing a prime factor with one of these special smoothness properties, and there are only "few" extensions of $\mathbb F_p$ small enough to give rise to efficient {$p-1$, $p+1$, etc.} methods like the above.

Amazingly, there's a brilliant workaround for this issue due to (Hendrik) Lenstra, known as the elliptic-curve method (ECM): Using elliptic curves instead of multiplicative groups of finite fields leaves us with exponentially more algebraic groups to choose from! One version of the basic strategy (with tons of possible optimizations — including the use of the Montgomery model, which was invented for this very purpose) proceeds as follows:

1. Pick two random elements $x,y$ of $\mathbb Z/n$.
2. Pick an elliptic curve $E$ "over $\mathbb Z/n$" such that $Q:=(x,y)$ is a point on $E$.
3. Pray to the deity of your choice that $Q$ is a point of smooth order modulo some prime divisor $p$ of $n$.
4. Compute $[e]Q$ for a large smooth integer $e$ as usual, simply acting as if $\mathbb Z/n$ was a field.
5. Now, if our prayers were heard, then the point $[e]Q$ actually turns into the point at infinity upon reduction modulo $p$. Hence, assuming the use of projective coordinates and writing $[x:y:z] = [e]Q$, simply computing $\gcd(z,n)$ at the end splits off $p$. (Recall that the point at infinity is $[0:1:0]$, the only point with $z$-coordinate $0$.)
   (For affine coordinates, what happens is that at some point we'll find ourselves trying to divide by an element $u$ of $\mathbb Z/n$ that's not coprime to $n$ — which means that $\gcd(u,n)$ splits $n$.)
6. If no divisor of $n$ was found, just try again (optionally switching to a different deity).

What have we gained? First, note that all the complexity of the algorithm lies in step 3: Smooth integers are rare, and we'll have to retry/pray a lot, but all the other steps of the algorithm are "fast". However, the order of $E$ modulo $p$ pretty much behaves like a uniformly random integer in the Hasse interval \[ [p+1-2\sqrt p ;\; p+1+2\sqrt p] \,\text, \] and the order of $Q$ won't be much smaller than the curve order. Assuming reasonably big $p$, the number of smooth integers in this interval depends only on the size of $p$, but (crucially) not on any other properties of $p$! Therefore, applying the strategy above gives rise to a general-purpose integer factorization algorithm whose complexity depends only on the size of the prime factors.

Having said all this, the potentially disappointing news is that the primes used in the f_cktoring challenge were chosen a little bit too large (≈ 256 bits) to be easily found by ECM within the 48-hour timeframe of the CTF. So, the question is:

How to hack?

With the discussion of algebraic-group-based factoring methods in mind, the solution to the challenge requires just one additional ingredient: In a nutshell, use ECM but choose wisely. That is, instead of trying random curves modulo $n$, we'll exploit our prior knowledge about $p$ to construct a curve that's guaranteed to have smooth order modulo $p$!

The key to choosing wisely is complex multiplication ("CM method"): Given two large integers $q,t$ such that the squarefree part of the discriminant \[ \Delta:=t^2-4q \] is tiny, one can efficiently construct an elliptic curve $E$ such that \[ \#E(\mathbb F_q) = q+1-t \,\text. \]

In our case, we'd like $q$ to equal the $p=(x{+}1)^2+123y^2$ from above and $\#E(\mathbb F_p)$ to be the smooth value $x^2+123y^2$. Thus \[ t = q + 1 - \#E(\mathbb F_p) = (x+1)^2+123y^2 + 1 - (x^2+123y^2) = 2x+2 \] and therefore \[ \Delta = t^2-4q = (2x+2)^2 - 4(x+1)^2 - 4\cdot 123y^2 % = 4x^2 + 8x + 4 - 4x^2 - 8x - 4 - 4\cdot 123y^2 = -123\cdot 4y^2 \,\text. \] Clearly $4y^2$ is a square, hence the squarefree part of $\Delta$ is simply $-123$. If that ain't small, I don't know what is!

Curves from $\Delta$

The theory underlying this is really beautiful and deep, but in practice, all we need is a Hilbert class polynomial: They are a family of univariate polynomials with integer coefficients, defined by the condition that they vanish precisely at the $j$-invariants of elliptic curves whose endomorphism ring is isomorphic to a particular imaginary quadratic ring. This may sound scary, but modulo $p$, the condition simply boils down to having a particular value of $\Delta$, and the size of $\Delta$ influences the complexity. Our $\Delta$ is huge, but the cardinality of the curve — which is all we need — actually depends only on the squarefree part of $\Delta$; thus, we can ignore the $4y^2$ and use $-123$ instead.

Concretely, the Hilbert class polynomial of the ring $\mathbb Z[\sqrt{-123}]$ is \[ H_{-123}(x) = x^2 + 1354146840576000 x + 148809594175488000000 \,\text, \] hence (by the calculation of $\Delta$ above) for any prime $p$ of the form $(x{+}1)^2+123y^2$, the elliptic curves over $\mathbb F_p$ defined by the roots of $H_{-123} \bmod p$ will have precisely $x^2+123y^2$ points. (Actually, this is a lie, but it's a minor technical detail: The cardinality might belong to the quadratic twist.)

Isn't this circular?

But y⁷, I hear you ask, how do we find a root of this polynomial modulo $n$ without factoring $n$ first? That's an excellent point.

Answer: We don't. Instead, compute everything with a symbolic solution to this equation. In other words, work in the quotient ring $\mathbb Z[x]/(n,H_{-123}(x))$, and use a curve with $j$-invariant $x$. The point is that any computations performed using the symbolic root $x$ must also hold for any "actual" root of $H_{-123}\bmod n$, and therefore everything will be okay in the end.

Code

Here's my solver (not optimized or beautified, but good enough):

#!/usr/bin/env sage
import re
proof.all(False)

d = 123
n,c = map(ZZ, re.findall('0x[0-9a-f]+', open('output.txt').read()))

# Sage can't invert in ℤ[x]/(n,f(x)); we use this workaround
def inv(f):
    ff = f.lift().change_ring(QQ)
    gg = f.parent().modulus().change_ring(QQ)
    hh = xgcd(ff,gg)            # compute xgcd over ℚ
    assert hh[0] == 1
    return f.parent()(hh[1])    # reduce back

# X1 is P-Q, X2 is P, X3 is Q, returns 2P, P+Q
def xDBLADD(a,b,X1,X2,X3):
    Z1 = Z2 = Z3 = 1
    if X1 == (): X1,Z1 = 1,0
    if X2 == (): X2,Z2 = 1,0
    if X3 == (): X3,Z3 = 1,0
    X4 = (X2^2-a*Z2^2)^2-8*b*X2*Z2^3
    Z4 = 4*(X2*Z2*(X2^2+a*Z2^2)+b*Z2^4)
    R = 2*(X2*Z3+X3*Z2)*(X2*X3+a*Z2*Z3)+4*b*Z2^2*Z3^2
    S = (X2*Z3-X3*Z2)^2
    X5 = R-S*X1
    Z5 = S
    return (X4*inv(Z4) if Z4 else ()), (X5*inv(Z5) if Z5 else ())

def xMUL(a,b,n,P):
    Q,PQ = (),P
    for bit in map(int,f'{n:b}'):
        if bit:
            P,Q = xDBLADD(a,b,PQ,P,Q)
        else:
            Q,P = xDBLADD(a,b,PQ,Q,P)
    return Q

################################################################

ls = list(sorted({u^2+d*v^2 for u in range(49+1) for v in range(7+1)} - {0}))
L = lcm(ls)     # or something

def fac(m, path=''):

    print(f'\x1b[33m    {"["+path+"]":9} ({int(m).bit_length()} bits)\x1b[0m')

    m = ZZ(m)
    if is_pseudoprime(m):                   # our job is done here
        print(f'--> \x1b[36m{m:#x}\x1b[0m'); print()
        return [m]

    R.<x> = Zmod(m)[]
    S.<J> = R.quotient(hilbert_class_polynomial(-d))                # J is a symbolic root of H_{-d} modulo n
    a0, b0 = (-3*J^2+3*1728*J), (-2*J^3+4*1728*J^2-2*1728^2*J)      # y^2 = x^3 + a0*x + b0 has j-invariant J

    for i in range(1,20):

        e = L^i                             # large smooth exponent
        P = J.parent().random_element()     # random point
        try:
            Q = xMUL(a0,b0,e, P)
            print(f'{i:2} nope')
#            print(P, Q)
            continue
        except ZeroDivisionError as e:      # probably found a divisor!
            vs = list(map(ZZ, re.findall('[0-9]+', e.args[0])))

        f = gcd(vs[0],m)
        if 1 < f < m:
            print(f'\x1b[34m{m:#x} = {f:#x} * {m//f:#x}\x1b[0m'); print()
            return fac(f, path+'L') + fac(m//f, path+'R')   # recursively factor the rest

    print('\x1b[31mFAILED\x1b[0m')
    exit(1)

################################################################

fs = fac(n)

print('    \x1b[1;4mfactorization\x1b[0m')
for f in fs:
    print(f'\x1b[36m{f:#x}\x1b[0m')
print()

privkey = int(~Mod(65537, prod(f-1 for f in fs)))
print(f'\x1b[35md = {privkey:#x}\x1b[0m'); print()

flag = ZZ(pow(c, privkey, n))
flag = bytes.fromhex(f'{flag:x}').decode(errors='replace')
print(f'--> \x1b[1;32m{flag}\x1b[0m'); print()

This script takes about a minute to factor $n$ and decrypt the flag:

    []        (1836 bits)
 1 nope
0xb96069105fd32b73363cdf8eae91b1f5ec165b1ffaa8b9b56af7414e8cbeac03b897e8f8bea3067d7cd2c4baa4d1cedff94c6918fb266039057843f1ee43bd6e3b0fcb8eed1e0945a4eacb4f4c720a6fffaaadb14b6ef83549868db5a0fb39f29b97a94e1449f4b56a7775fc10cbe9e2a0d39c74ffcde60bf66f8b880654d4eb526283ab37a72e03fe1a68af0d2017ae04f7d3be6c2fad768d76a67ef2ceacac2b846ec3393f6c962861d06105134c5b7884bf561401a9fd897fff7aeb4c591a3cbd9f72b68eb4ae0ac64f213ef38449ca7f6ae98113a95731525661d5bab460b1b9f3eb801 = 0x1d3be7d70e0871d6df6533c9b811c79bc4b66347bfb7900e8a2b9303f070fa801 * 0x657543999adfed3d76740b8ef51ca29305585d7cd1ceae051be6fc8a06ed3ea9971e42a11c409a4fed38798ae400be86f1326d5f8aa731d4c5b2802b94f151a3ec09e84eb365c2edae2e095a1e8926f5f040556c55f2b430cdb2b6993687f0382af6961291fe9fa1935e953584a20b8cf86758666672f5f8683b33f586dc4c781460e4bb2f9e68ffa2ed0d181576c63283847831b5c2a03ca6d13eedc05cb00d1ee651afcaab34e3757b9ab8b44f94ddefc1b58a8430e03e93a037926264b6a326fc5af1001

    [L]       (257 bits)
--> 0x1d3be7d70e0871d6df6533c9b811c79bc4b66347bfb7900e8a2b9303f070fa801

    [R]       (1579 bits)
 1 nope
0x657543999adfed3d76740b8ef51ca29305585d7cd1ceae051be6fc8a06ed3ea9971e42a11c409a4fed38798ae400be86f1326d5f8aa731d4c5b2802b94f151a3ec09e84eb365c2edae2e095a1e8926f5f040556c55f2b430cdb2b6993687f0382af6961291fe9fa1935e953584a20b8cf86758666672f5f8683b33f586dc4c781460e4bb2f9e68ffa2ed0d181576c63283847831b5c2a03ca6d13eedc05cb00d1ee651afcaab34e3757b9ab8b44f94ddefc1b58a8430e03e93a037926264b6a326fc5af1001 = 0x5963a5406f0ba2cc9c2fbff49229511a631ef487ecc32433cdb2a417347fee80001 * 0x122904b76c14ac77a4d906dae6cd0e88d379933e6f2140c23e3d6e835a976a76a26a862c0a3ed0e45ac59285992b21b76691ca5bf5f9e8a85df90ddbf3bbeda73ddc2945c00aab1ac7530c500ca522dc39af4e9b74ae710607d7d489cbaa12658d266c7792718105c050a521e2f8673db9ddec4299faf8c8f9ca1ff44c75a74159155daa17d99bbce8451c32ad1c7a888dd43988e6bdd9abb605d18c97361a8e146c71001

    [RL]      (267 bits)
--> 0x5963a5406f0ba2cc9c2fbff49229511a631ef487ecc32433cdb2a417347fee80001

    [RR]      (1313 bits)
 1 nope
0x122904b76c14ac77a4d906dae6cd0e88d379933e6f2140c23e3d6e835a976a76a26a862c0a3ed0e45ac59285992b21b76691ca5bf5f9e8a85df90ddbf3bbeda73ddc2945c00aab1ac7530c500ca522dc39af4e9b74ae710607d7d489cbaa12658d266c7792718105c050a521e2f8673db9ddec4299faf8c8f9ca1ff44c75a74159155daa17d99bbce8451c32ad1c7a888dd43988e6bdd9abb605d18c97361a8e146c71001 = 0x2163cbb35018f93c31ce495921ef2ca081a9e8ce6f2b5fb5a6691cf47499ae001 * 0x8b3c0d3797a5aa61ae9ba1cfac61dc19c4cee522ce67976218a367561e6afdf9c25693bacf41b67b844d58a599288e5388a1d591848332d54bd0ce3b95f2b5bbdf52ba673d908e107ef8dd1fad747c134372495d81f1b6fab4d2a0648f94c31736a595aba50f29c2abc77aa70ca3063c4799f056080f9ec357a526e525f9028e732c3001

    [RRL]     (258 bits)
--> 0x2163cbb35018f93c31ce495921ef2ca081a9e8ce6f2b5fb5a6691cf47499ae001

    [RRR]     (1056 bits)
 1 nope
0x8b3c0d3797a5aa61ae9ba1cfac61dc19c4cee522ce67976218a367561e6afdf9c25693bacf41b67b844d58a599288e5388a1d591848332d54bd0ce3b95f2b5bbdf52ba673d908e107ef8dd1fad747c134372495d81f1b6fab4d2a0648f94c31736a595aba50f29c2abc77aa70ca3063c4799f056080f9ec357a526e525f9028e732c3001 = 0x509a9134d45e8a12743f619c7dfef456d3d619a2c60c6838629febe2b6c8ce40001 * 0x1ba3691cbe7eb58ca2593647558eb5f7ec5c3022b795bf3e25f3191f728494a32d410fbd4eb50a90f71d572ca53742fcdb60be2bec77a4a5e2237ca13e5b2a76160dc21ca78e0c32985b29cc0766cd901ca8ddd0a768d880371d4a8129049726483001

    [RRRL]    (267 bits)
--> 0x509a9134d45e8a12743f619c7dfef456d3d619a2c60c6838629febe2b6c8ce40001

    [RRRR]    (789 bits)
 1 nope
 2 nope
0x1ba3691cbe7eb58ca2593647558eb5f7ec5c3022b795bf3e25f3191f728494a32d410fbd4eb50a90f71d572ca53742fcdb60be2bec77a4a5e2237ca13e5b2a76160dc21ca78e0c32985b29cc0766cd901ca8ddd0a768d880371d4a8129049726483001 = 0x16fa7fa06d047d91c4ad338fa2eaf8f8a289b5d8655d937bf2f2512f82a083001 * 0x133ea0c20f32625834176be9fa86d1b2e4a1513257beb84639efee0f2e68a096127c25d73f0f7409e2ce7a371cd886ee286a2e382a622de408488df3fa0052fc400001

    [RRRRL]   (257 bits)
--> 0x16fa7fa06d047d91c4ad338fa2eaf8f8a289b5d8655d937bf2f2512f82a083001

    [RRRRR]   (533 bits)
 1 nope
0x133ea0c20f32625834176be9fa86d1b2e4a1513257beb84639efee0f2e68a096127c25d73f0f7409e2ce7a371cd886ee286a2e382a622de408488df3fa0052fc400001 = 0x1a9568e26c70252463c0e795e08c053220b35924678285ecef785ed67d055e00001 * 0xb9535bb06458d9a6995594bc4077f751063cabd7f1ca47edcff4c7d2482a6600001

    [RRRRRL]  (265 bits)
--> 0x1a9568e26c70252463c0e795e08c053220b35924678285ecef785ed67d055e00001

    [RRRRRR]  (268 bits)
--> 0xb9535bb06458d9a6995594bc4077f751063cabd7f1ca47edcff4c7d2482a6600001

    factorization
0x1d3be7d70e0871d6df6533c9b811c79bc4b66347bfb7900e8a2b9303f070fa801
0x5963a5406f0ba2cc9c2fbff49229511a631ef487ecc32433cdb2a417347fee80001
0x2163cbb35018f93c31ce495921ef2ca081a9e8ce6f2b5fb5a6691cf47499ae001
0x509a9134d45e8a12743f619c7dfef456d3d619a2c60c6838629febe2b6c8ce40001
0x16fa7fa06d047d91c4ad338fa2eaf8f8a289b5d8655d937bf2f2512f82a083001
0x1a9568e26c70252463c0e795e08c053220b35924678285ecef785ed67d055e00001
0xb9535bb06458d9a6995594bc4077f751063cabd7f1ca47edcff4c7d2482a6600001

d = 0x8aaeb2cb1a353beb4c674c044d5b181cd7cb7dceb5d266646faa9c83ec4d7517c00fc47f755d07b7261a222083f0d52f209559b2901c1acca38a34cfa06b4a3aa1e848b3b74f7ecccc78769c4826f3ad10684e158a6cf0d795ec16986a62ace7a81097bd7b51433b9b50fd3ebf7db1342ea4a1486f6aef8b937fa180c8bb184836895e19043f5fa39758e8cf3074f040b7bd91dfebfbee10910e354f16cb4d6da8eff222810a24331cbe1625a50712f92ebdd0cbdc6686ab3df3167fbcc9cae2ee053724353a8bb3018dd6fbfa1babb2cfaa46bdb4b99670000ffff0000ffff0000ffff0001

--> hxp{M4tH_1s_l1T3r4Lly_M4GiC}

A twist of fate

The challenge was made in July, and I hadn't found this technique in the literature anywhere. I was living more or less happily as the months were passing, for little did I know that A WEEK BEFORE THE CTF, this paper by Giuseppe Vitto would appear online, which not only explains in meticulous detail how to solve the problem, but even goes as far as providing an implementation which I've been told breaks the challenge with only minor changes. Worst timing ever, and it's all because Giuseppe recklessly skipped submitting the paper to hxp's responsible-disclosure process first!

I hope people still had fun thinking about this problem, despite the fact that it kind of degraded to an implement-the-paper or even just run-the-code-properly task for people who found the paper. Luckily(?), the solve count suggests that the set of such people is rather small.

In addition, it turns out that the same method had been published at least two times before, first in 2002 by Cheng and then in 2019 by Aikawa, Nuida, and Shirase. Given the connections to $p\pm1$ and ECM, I think it's not too surprising that this idea occured to multiple people in the history of mankind.

Moral of the story: idk 🤷‍♂️