solve script for Curta Challenge 8, "Groovy"

solve.sage

#!/usr/bin/env sage
punch = 906459278810089239293436146013992401709

# challenge value here
s = 13763640145752339203207

def verify(s, a,b):
    """solution verifier ported from the contract"""
    c = 1
    assert a > 100 and b > 100
    assert smoothie(a,b) == 1
    assert (a * milkshake(b + c) % punch + b * milkshake(a + c) % punch + c * milkshake(a + b) % punch) % punch == s

def blender(a,b):
    """turns out this is the extended GCD"""
    if b == 0:
        return a,1, 0
    y, aj, bj = blender(b, a%b)
    return y, bj, aj - a//b * bj

def milkshake(a):
    """this is the modular multiplicative inverse"""
    _,aj,_ = blender(a, punch)
    assert aj > 0 #lol CTS is too lazy to handle all cases
    assert a*aj % punch == 1
    return aj % punch

def smoothie(a, b):
    """and this is the gcd"""
    y,_,_ = blender(a,b)
    return y

"""looking at milkshake, we are working a prime field"""
GF = GF(punch,1)


"""
We have two unknowns with a cubic relation between them.
We use a rejection sampling approach where we take one at random
and solve for the other one.
"""
while True:
    R.<a> = GF[]
    b = randint(100, 0x7fffffffffffffffffffffffffffffff)
    c = 1
    cubic = (a/(b+c) + b/(a+c) + c/(a+b) - s).numerator()
    roots = cubic.roots()
    if not roots:
        continue
    a = int(roots[0][0])
    if not 100 <= a <= 0x7fffffffffffffffffffffffffffffff:
        continue
    if math.gcd(int(a), int(b)) != 1:
        continue
    try:
        verify(s, a,b)
    except AssertionError:
        continue
    break

print(f'{a=}')
print(f'{b=}')
solution = int(b) << 128 | int(a)
print(hex(solution))