1
# !/usr/bin/env sage
2
import re
3
from math import gcd
4
from gmpy2 import iroot
5

6
def parse_output(path):
7
    vals = {}
8
    points = {}
9
    with open(path, "r", encoding="utf-8") as f:
10
        for line in f:
11
            line = line.strip()
12
            if not line or "=" not in line:
13
                continue
14
            k, v = [x.strip() for x in line.split("=", 1)]
15
            if k[0] in ("x", "y", "r") and k[1:].isdigit():
16
                idx = int(k[1:])
17
                points.setdefault(idx, {})[k[0]] = Integer(v)
18
            else:
19
                vals[k] = Integer(v)
20
    return vals, points
21

22
def interpolate_mod_n(points, n):
23
    Zn = Zmod(n)
24
    R = PolynomialRing(Zn, "X")
25
    X = R.gen()
26

27
    f = R(0)
28
    m = len(points)
29
    for i in range(m):
30
        xi, yi = points[i]
31
        num = R(1)
32
        den = Zn(1)
33
        for j in range(m):
34
            if i == j:
35
                continue
36
            xj, _ = points[j]
37
            num *= (X - xj)
38
            den *= (xi - xj)
39

40
        den_int = int(den)
41
        g = gcd(den_int, n)
42
        if g != 1:
43
            raise ValueError(f"Non-invertible interpolation denominator; gcd={g}")
44
        f += yi * num * den.inverse_of_unit()
45
    return f
46

47
def gcd_over_zn(a, b, n):
48
    R = a.parent()
49
    while b != 0:
50
        lc = b.leading_coefficient()
51
        g = gcd(int(lc), n)
52
        if g != 1:
53
            return None, g
54
        b = b * lc.inverse_of_unit()
55
        a, b = b, a % b
56

57
    if a == 0:
58
        return R(0), None
59

60
    lc = a.leading_coefficient()
61
    g = gcd(int(lc), n)
62
    if g != 1:
63
        return None, g
64
    a = a * lc.inverse_of_unit()
65
    return a, None
66

67
def int_to_bytes(x):
68
    x = int(x)
69
    if x == 0:
70
        return b"\x00"
71
    return x.to_bytes((x.bit_length() + 7) // 8, "big")
72

73
def main():
74
    vals, raw_points = parse_output("output.txt")
75
    n = int(vals["n"])
76
    e = int(vals["e"])
77
    c1 = int(vals["c1"])
78
    c2 = int(vals["c2"])
79

80
    points = []
81
    for idx in sorted(raw_points):
82
        x = int(raw_points[idx]["x"])
83
        y_noisy = int(raw_points[idx]["y"])
84
        r_enc = int(raw_points[idx]["r"])
85

86
        r, exact = iroot(r_enc, e)
87
        if not exact:
88
            raise ValueError(f"r{idx} root is not exact")
89
        r = int(r)
90
        y = (y_noisy - r) % n
91
        points.append((x, y))
92

93
    f = interpolate_mod_n([(Zmod(n)(x), Zmod(n)(y)) for x, y in points], n)
94

95
    R = f.parent()
96
    X = R.gen()
97
    P = X**e - Zmod(n)(c1)
98
    Q = f**e - Zmod(n)(c2)
99

100
    g, nontrivial = gcd_over_zn(P, Q, n)
101

102
    m = None
103
    if g is not None and g.degree() == 1:
104
        b = g[0]
105
        m = int((-b) % n)
106
    elif nontrivial is not None and 1 < nontrivial < n:
107
        p = int(nontrivial)
108
        q = n // p
109
        phi = (p - 1) * (q - 1)
110
        d = inverse_mod(e, phi)
111
        m = int(pow(c1, int(d), n))
112
    else:
113
        gg = P.gcd(Q)
114
        if gg.degree() == 1:
115
            m = int((-gg[0] / gg[1]) % n)
116

117
    if m is None:
118
        raise RuntimeError("Failed to recover m")
119

120
    msg = int_to_bytes(m)
121
    print("m_bytes:", msg)
122

123
    try:
124
        s = msg.decode()
125
    except Exception:
126
        s = msg.decode(errors="ignore")
127

128
    print("decoded:", s)
129
    mm = re.search(r"CodeVinci\{[^}]+\}", s)
130
    if mm:
131
        print("FLAG:", mm.group(0))
132

133
if __name__ == "__main__":
134
    main()