# Exact checks (Gaussian rationals) of the SL2/SO2 computations: stabilizer of x_m (A2-1),
# the slice at x_1 inside closure(O_2) and the slice at x_0 inside closure(O_2) (A2-6).
from fractions import Fraction as Fr
import random
random.seed(7)

class G:  # Gaussian rational a + b i
    def __init__(s, a, b=0): s.a, s.b = Fr(a), Fr(b)
    def __add__(s, o): o = o if isinstance(o, G) else G(o); return G(s.a + o.a, s.b + o.b)
    __radd__ = __add__
    def __neg__(s): return G(-s.a, -s.b)
    def __sub__(s, o): return s + (-(o if isinstance(o, G) else G(o)))
    def __rsub__(s, o): return G(o) - s
    def __mul__(s, o): o = o if isinstance(o, G) else G(o); return G(s.a*o.a - s.b*o.b, s.a*o.b + s.b*o.a)
    __rmul__ = __mul__
    def __eq__(s, o): o = o if isinstance(o, G) else G(o); return s.a == o.a and s.b == o.b
    def iszero(s): return s.a == 0 and s.b == 0
I = G(0, 1)
def rnd(): return G(Fr(random.randint(-9, 9), random.randint(1, 5)), Fr(random.randint(-9, 9), random.randint(1, 5)))

# Laurent polynomials in t: dict exponent -> G
def lmul(x, y):
    z = {}
    for i, u in x.items():
        for j, v in y.items():
            z[i + j] = z.get(i + j, G(0)) + u * v
    return {k: v for k, v in z.items() if not v.iszero()}
def ladd(x, y, s=1):
    z = dict(x)
    for k, v in y.items(): z[k] = z.get(k, G(0)) + (v if s == 1 else -v)
    return {k: v for k, v in z.items() if not v.iszero()}
def mm(A, B): return [[ladd(lmul(A[i][0], B[0][j]), lmul(A[i][1], B[1][j])) for j in range(2)] for i in range(2)]
def det(A): return ladd(lmul(A[0][0], A[1][1]), lmul(A[0][1], A[1][0]), -1)
def tr(A): return [[A[0][0], A[1][0]], [A[0][1], A[1][1]]]
def minval(A): return min(min(e.keys()) for row in A for e in row if e)

ok = True
# (1) stabilizer: g = [[a, b t^m],[-b t^-m, a]] with a^2+b^2 = 1 fixes x_m = diag(t^m, t^-m) exactly
for m in (1, 2, 3):
    # rational point on a^2+b^2=1: a=(1-s^2)/(1+s^2), b=2s/(1+s^2), times t^m for b (power series not needed: exact identity)
    s = rnd()
    den = G(1) + s * s
    if den.iszero(): continue
    # use a, b with a^2+b^2=1 as Gaussian rationals (constant in t): identity is polynomial in a, b
    from fractions import Fraction
    # invert den
    n2 = den.a * den.a + den.b * den.b
    inv = G(den.a / n2, -den.b / n2)
    a = (G(1) - s * s) * inv; b = G(2) * s * inv
    assert (a * a + b * b) == G(1)
    g = [[{0: a}, {m: b}], [{-m: -b}, {0: a}]]
    xm = [[{m: G(1)}, {}], [{}, {-m: G(1)}]]
    ok &= mm(mm(g, xm), tr(g)) == [[{m: G(1)}, {}], [{}, {-m: G(1)}]] and det(g) == {0: G(1)}
print('stabilizer identity g x_m g^T = x_m, det g = 1:', ok)

# (2) slice at x_1: h = [[1 + c/t + c^2/t^2 + c^3/t^3, e/t],[e/t^3, 1 - c/t]], e^2 = -c^4
ok2 = True
for _ in range(20):
    c = rnd()
    for e in (I * c * c, -I * c * c):
        h = [[{0: G(1), -1: c, -2: c * c, -3: c * c * c}, {-1: e}], [{-3: e}, {0: G(1), -1: -c}]]
        x1 = [[{1: G(1)}, {}], [{}, {-1: G(1)}]]
        A = mm(h, x1)
        sym = A[0][1] == A[1][0]
        ok2 &= det(h) == {0: G(1)} and sym and minval(A) == -2
        # h in L^{<0}G and t^{-lam} h t^{lam} in L^{<0}G (entries: diag 1 + t^{-1}C[t^{-1}], offdiag t^{-1}C[t^{-1}])
        conj = [[h[0][0], {k - 2: v for k, v in h[0][1].items()}], [{k + 2: v for k, v in h[1][0].items()}, h[1][1]]]
        for M in (h, conj):
            for i in range(2):
                for j in range(2):
                    for k, v in M[i][j].items():
                        ok2 &= (k < 0) or (k == 0 and i == j and v == G(1))
print('slice at x_1: two branches e = +-i c^2 lie in closure(O_2), meet only at x_1:', ok2)

# general element of that slice: h = [[1+a1 u+a2 u^2+a3 u^3, b u],[b u^3, 1+d u]], u = 1/t; det = 1 forces
# a1=-d, a2=d^2, a3=-d^3, b^2=-d^4 (coefficients of u, u^2, u^3, u^4).  Check the coefficient formulas:
def detcoeffs(a1, a2, a3, b, d):
    # (1+a1u+a2u^2+a3u^3)(1+du) - b^2 u^4
    return [a1 + d, a2 + a1 * d, a3 + a2 * d, a3 * d - b * b]
ok3 = True
for _ in range(20):
    d = rnd(); b = I * d * d
    ok3 &= all(x.iszero() for x in detcoeffs(-d, d * d, -d * d * d, b, d))
print('determinant conditions at x_1 solved by the two-branch family:', ok3)

# (3) slice at x_0: h = [[1+a1 u+a2 u^2, b1 u + b2 u^2],[same, 1 - a1 u + d2 u^2]]: equations
def eqs(a1, b1, a2, b2, d2):
    return [a2 + d2 - a1 * a1 - b1 * b1, a1 * (d2 - a2) - G(2) * b1 * b2, a2 * d2 - b2 * b2]
def deth(a1, b1, a2, b2, d2):
    h = [[{0: G(1), -1: a1, -2: a2}, {-1: b1, -2: b2}], [{-1: b1, -2: b2}, {0: G(1), -1: -a1, -2: d2}]]
    return det(h)
ok4 = True
for _ in range(20):
    P = rnd(); Q = rnd()
    if P.iszero(): continue
    n2 = P.a * P.a + P.b * P.b; Pinv = G(P.a / n2, -P.b / n2)
    R = Q * Q * Pinv
    pt = (G(2) * Q, R - P, (P + R) * P, (P + R) * Q, (P + R) * R)   # Y_main via (P,Q,R), PR = Q^2
    ok4 &= all(x.iszero() for x in eqs(*pt)) and deth(*pt) == {0: G(1)}
    a1 = rnd(); nu = rnd()
    for sg in (1, -1):
        pt = (a1, sg * I * a1, nu, sg * I * nu, -nu)                   # planes P_+, P_-
        ok4 &= all(x.iszero() for x in eqs(*pt)) and deth(*pt) == {0: G(1)}
print('slice at x_0: cone {PR=Q^2} and planes P_+, P_- satisfy det h = 1:', ok4)
# intersections: Y_main meets P_+ along the line (a1, i a1, 0,0,0) (points with P+R = 0, i.e. R = -P, Q = +-iP)
P = rnd(); Q = I * P; R = -P
pt = (G(2) * Q, R - P, (P + R) * P, (P + R) * Q, (P + R) * R)
print('Y_main point with R=-P, Q=iP:', [ (x.a, x.b) for x in pt ], ' -> b1 = i a1 ?', pt[1] == I * pt[0])