hxp 38C3 CTF 2024 - alcoholic_variety (crypto)
Published on 30 Dec 2024
Writeups
This challenge was about computing discrete logarithm on an obfuscated variety given by a 9-dimensional projective embedding.
Server code:
sec = SystemRandom().randrange(p**2)
pub = scalarmul(base, sec)
print(f'{pub = }')
sol = int(input(), 0)
if not 0 <= sol < p**2:
exit(1)
chk = scalarmul(base, sol)
if chk != pub:
exit(1)
Variety arithmetic:
def add(P, Q):
X0, X1, X2, X3, X4, X5, X6, X7, X8, X9 = P
Y0, Y1, Y2, Y3, Y4, Y5, Y6, Y7, Y8, Y9 = Q
T00 = X0 * Y0 % p
T01 = X0 * Y1 % p
...
S0000 = T00 * T00 % p
S0001 = T00 * T01 % p
...
Z0 = 0x5eb3705ee948c25660e2bf8781a1dbf6c54a6bdb * S0000
Z0 += 0x3fed1b949b315aef4ef4021c37edba708b9a6919 * S0001
Z0 += 0x8a49fca770fdcc94c153d9c95cada504c13ca728 * S0002
...
Z1 = 0xbdfaabf8789d6aff2614fecf066d394d07b195ea * S0000
Z1 += 0x4492da5c2027c6ee0b16e746b34fa2254278a429 * S0001
Z1 += 0x428294f92b85f159d25672db96204ccf21db67e2 * S0002
...
return Z0, Z1, Z2, Z3, Z4, Z5, Z6, Z7, Z8, Z9
Interpolate projection on 2 affine coordinates
Let's try to see if any two coordinates of the variety are constrained by a low degree polynomial. This should distinguish a curve (1-dimensional variety) from higher-dimensional objects. We can do this by naive interpolation.
In [9]:
from vuln import *
In [10]:
deg = 9 # total degree
monos = []
for i in range(deg + 1):
for j in range(deg + 1 - i):
monos.append((i, j))
gen = base
pt = gen
data = []
for i in range(2, len(monos)+5):
pt = normalize(add(pt, gen))
# first 2 affine coordinates (since normalized)
xs = pt[1:3]
row = [
prod(x**e for x, e in zip(xs, es))
for es in monos
]
data.append(row)
mat = matrix(GF(p), data)
rk = mat.right_kernel()
rk
Out [10]:
Vector space of degree 55 and dimension 1 over Finite Field of size 1349317116055496347935183323829721306283718936287 Basis matrix: 1 x 55 dense matrix over Finite Field of size 1349317116055496347935183323829721306283718936287
We have 1 relation (curve)! Let's reconstruct the polynomial:
In [11]:
xs = PolynomialRing(GF(p), names="x, y").gens()
for sol in rk.matrix():
assert len(sol) == len(monos)
poly = []
for c, es in zip(sol, monos):
mono = prod(x**e for x, e in zip(xs, es))
poly.append(c * mono)
poly = sum(poly)
print(poly)
1290422967115358566266590777853979389598846205666*x^9 + 577379728613493404559030761507813969909021892646*x^8*y + 1243591010103894033874931036346746932772777061809*x^7*y^2 + 802820083981430427359276555319159610338384321310*x^6*y^3 + 1232868555484873001243395461503072086003412440799*x^5*y^4 + 1137483247688794514748538208782775965405225990619*x^4*y^5 + 1098058599413970905114088603603094763519900428362*x^3*y^6 + 521491890352857306282988928668018529753246687918*x^2*y^7 + 1260410534149648293496028865098499989940595729366*x*y^8 + 842782215788949538115688200389919742540609803787*y^9 + 148652987333517774706626469674067058579773949330*x^8 + 704072252718806103236352463412377087321883390307*x^7*y + 969550413402677076033559674669970647059658559811*x^6*y^2 + 914504036011531344190663463577438459411937411774*x^5*y^3 + 1033218420556944139578595749588116845242280409381*x^4*y^4 + 184031932531709285127717045341600269148663709973*x^3*y^5 + 522977194194553309971023227001628408398393067665*x^2*y^6 + 526942879502285241403569507339108847069792634647*x*y^7 + 17009683443167968083917508140397063144005599766*y^8 + 1206984109122888717952309707979157801919379180440*x^7 + 666646221516501167407644038840016851578513381070*x^6*y + 951883726405392312011964097975105374209885276989*x^5*y^2 + 645823893885258859236272557206945939990660011679*x^4*y^3 + 652478917929683952185569824651099051781512339594*x^3*y^4 + 1147786275280675470977066928226574352664312025618*x^2*y^5 + 511780223753123772491966598541792472900715927890*x*y^6 + 1242714503544177641576865227346261839104821592387*y^7 + 1119608915578092664931027342305520214618342542817*x^6 + 219465047216145355211916063856654026557709438956*x^5*y + 440868630139328723390887666440338107454072029138*x^4*y^2 + 884841851333430383561092868693582396236627016059*x^3*y^3 + 1074014609963442522186489331612835495173528091285*x^2*y^4 + 632915613590123123726310817491414055711594675603*x*y^5 + 360153763588701435270187378411040859528860441053*y^6 + 1312967192351481772025312036859334531188871704480*x^5 + 370542169622946669850165726450244003690345761279*x^4*y + 656709532575879601531214822345460775013507981650*x^3*y^2 + 511896697232660276092068265004006902231292911556*x^2*y^3 + 1070096180258818183515005064468343611321309367240*x*y^4 + 263411577279111351351767233648575110516269742212*y^5 + 713108719042207102637623462493952328733444108683*x^4 + 1060825606481194737492064176412987247144297085300*x^3*y + 822052194696345751139804593504655648005104702647*x^2*y^2 + 1114929232706754066135276260250399322468687032746*x*y^3 + 195141015153221254035344067042632562110836626787*y^4 + 1015616695064203722876228358981056451279442812366*x^3 + 282436498401508215802299272916904057251062041895*x^2*y + 271013532598523925647618619209267801341794393684*x*y^2 + 1311804110854716862303352250687716364825835352313*y^3 + 313489185139358333253607840927329962919747187302*x^2 + 30541526957092351611848006450395151067734672851*x*y + 1011318165555535557007129300549155976507803763648*y^2 + 379343306669746198394996028074861013023099496937*x + 902827703943354025101323492246819516991375781791*y + 1
Note: the method above works for any pair of coordinates, yielding the same result (up to isomorphism).
Magma magic
Let's check the curve with Magma. What genus does it have?
In [13]:
INTRO = f"""
p := 0xec5978da7a46c1326123bddc340f0b51710706df;
Fp := GF(p);
P<x,y,h> := ProjectiveSpace(GF(p), 2);
C := Curve(P, {poly.homogenize()});
"""
magma_free(INTRO + "Genus(C);")
Out [13]:
1
It's just an elliptic curve! We can use this elliptic curve constructor in Magma.
In [19]:
def parseEll(line):
line = line.split(" by ")[1].split(" over ")[0]
R = PolynomialRing(GF(p), names='x,y')
line = line.replace("=", "-(") + ")"
pol= R(line)
return EllipticCurve(pol)
zero = normalize(zero)
INTRO_ELL = INTRO + f"""
zero := C![Fp!{zero[1]},Fp!{zero[2]},Fp!1];
E, morph := EllipticCurve(C, zero);
DP := DefiningPolynomials(morph);
"""
R = PolynomialRing(GF(p), names='x,y,h')
# due to output size limitation,
# recover maps one per request
E = fx = fy = fz = None
for what in ("e", "fx", "fy", "fz"):
print(what)
out = dict(
e="E",
fx="DP[1]",
fy="DP[2]",
fz="DP[3]",
)[what]
res = magma_free(INTRO_ELL + out)
print(" ", res[:50], " ... ", res[-50:])
if what == "e":
E = parseEll(res.strip())
elif what == "fx":
fx = R(res.strip())
elif what == "fy":
fy = R(res.strip())
elif what == "fz":
fz = R(res.strip())
e
Elliptic Curve defined by y^2 + 134845292611396518 ... 1349317116055496347935183323829721306283718936287)
fx
369185259345838176777052263906592997715549392180*x ... 395550858285139424016691192951431539250724931*h^30
fy
1212997578007365840857289354447777066829450093798* ... 853383776805192036318629666051564276006794873*h^30
fz
y^30 + 8795495203214332856229038197518207798026256 ... 036696400080331365819900914614310084861522230*h^30
Elliptic Curve
Let's check the curve's order:
In [20]:
E
Out [20]:
Elliptic Curve defined by y^2 + 1348452926113965188163925155490560842831053962286*x*y + 322542631082841025630294763070556383220924470645*y = x^3 + 246451750546204524575955295721790114141002989457*x^2 + 507792243786200675949807212210783619593791522364*x + 599607646593483263696061470594911408845255474935 over Finite Field of size 1349317116055496347935183323829721306283718936287
In [22]:
factor(E.order())
Out [22]:
31306060063 * 6565122100096804799^2
Let's check the base point order:
In [44]:
def toE(pt):
return E(fx(*pt[1:3], 1), fy(*pt[1:3], 1), fz(*pt[1:3], 1))
G = toE(base)
factor(G.order())
Out [44]:
31306060063 * 6565122100096804799
In [26]:
scalarmul(base, 31306060063 * 6565122100096804799) == zero
Out [26]:
True
This means that the base point on the variety has the same order as the base point on the elliptic curve. Therefore, we can ignore the other coordinates from the AV.
ECDLP
The subgroup of size 31306060063 is small enough for BSGS, while for 6565122100096804799 it is easy to check that MOV is applicable, since full 6565122100096804799-torsion is present in GF(p).
In [30]:
n1 = 31306060063
n2 = 6565122100096804799
while True:
G0 = G*n1
GG = E.random_element()*n1
if GG.order() == n2:
if G0.weil_pairing(GG, n2) != 1:
break
def weil(R):
return (R*n1).weil_pairing(GG, n2)
wG = weil(G)
e = randrange(p**2)
assert weil(G)**e == weil(toE(scalarmul(base, e)))
wG
Out [30]:
59778686898639018472079341752054945291009307440
In [31]:
p
Out [31]:
1349317116055496347935183323829721306283718936287
$ ./cado-nfs.py --dlp --ell 6565122100096804799 1349317116055496347935183323829721306283718936287
...
Info:root: If you want to compute one or several new target(s), run ./cado-nfs.py /tmp/cado.sa7vsbno/p50.parameters_snapshot.0 target=<target>[,<target>,...]
$ ./cado-nfs.py /tmp/cado.sa7vsbno/p50.parameters_snapshot.0 target=59778686898639018472079341752054945291009307440
...
5226351864318504434
In [32]:
logG = 5226351864318504434
$ nc 195.201.94.102 17017
pub = (1, 735406580662044781122885749242460238340586667300, 941722307540000824428771210612263684979121415575, 527326145907634170351473774050826602536246641349, 825310686595543556039897141091632508776643285955, 1154484822650118998943388626264957946707287231804, 495656494315470026517165259848292961486617526527, 400252245628339165498755320503996683996311988179, 640462851595800008949362550256757141698758100550, 784884675235802092004177258420482016420009917756)
In [50]:
pub = (1, 735406580662044781122885749242460238340586667300, 941722307540000824428771210612263684979121415575, 527326145907634170351473774050826602536246641349, 825310686595543556039897141091632508776643285955, 1154484822650118998943388626264957946707287231804, 495656494315470026517165259848292961486617526527, 400252245628339165498755320503996683996311988179, 640462851595800008949362550256757141698758100550, 784884675235802092004177258420482016420009917756)
weil(toE(pub))
Out [50]:
74395282880777674911447620422506821657920403862
$ ./cado-nfs.py /tmp/cado.sa7vsbno/p50.parameters_snapshot.0 target=74395282880777674911447620422506821657920403862
...
3813075515379119385
In [54]:
logQ = 3813075515379119385
Finishing ECDLP
BSGS part (using Sage):
In [53]:
Q = toE(pub)
sec1 = (Q * n2).log(G * n2)
sec1 = int(sec1) % n1
print("sec1", sec1)
sec1 28150040067
MOV part:
In [57]:
sec2 = GF(n2)(logQ) / GF(n2)(logG)
sec2
Out [57]:
4592982170174041327
In [59]:
ans = crt([int(sec1), int(sec2)], [n1, n2])
print(scalarmul(base, ans) == pub)
ans
True
Out [59]:
198281955505251745732602763692
$ nc -vn 195.201.94.102 17017 (base)
(UNKNOWN) [195.201.94.102] 17017 (?) open
pub = (1, 735406580662044781122885749242460238340586667300, 941722307540000824428771210612263684979121415575, 527326145907634170351473774050826602536246641349, 825310686595543556039897141091632508776643285955, 1154484822650118998943388626264957946707287231804, 495656494315470026517165259848292961486617526527, 400252245628339165498755320503996683996311988179, 640462851595800008949362550256757141698758100550, 784884675235802092004177258420482016420009917756)
(continued)
198281955505251745732602763692
hxp{Wh4t_1s_c0mmUt4T1v3_4Nd_h1gH_d1M3Ns10n4L?}