"""Directed interval controllability checks adapted from the paper."""
from mpmath import mp,iv
from local_certificate import Jet,source

def jac(radius,cert):
    box = [iv.mpf(x) + iv.mpf(['-' + radius, radius]) for x in cert['root_center']]
    jets, _ = source(Jet(box[0]), Jet(box[1]))
    return iv.matrix([[x.v.real for x in row] for row in jets])


def kalman(J, b):
    v = iv.matrix(b); cols = []
    for k in range(9):
        cols.append(v); v = J*v/100
    A = [[cols[j][i] for j in range(9)] for i in range(9)]
    piv, sign = [], 1
    for k in range(9):
        p = max(range(k, 9), key=lambda i: float(abs(A[i][k]).mid))
        if p != k:
            A[k], A[p] = A[p], A[k]; sign = -sign
        pv = A[k][k]
        if not (pv.a > 0 or pv.b < 0):
            return None, piv
        piv.append(pv)
        for i in range(k + 1, 9):
            mlt = A[i][k]/pv
            for j in range(k + 1, 9):
                A[i][j] -= mlt*A[k][j]
            A[i][k] = iv.mpf(0)
    det = iv.mpf(sign)
    for pv in piv:
        det *= pv
    return det, piv


def e(*idx_sign):
    v = [0]*9
    for i, s in idx_sign:
        v[i] = s
    return v


DIRS = {}
for i in range(3):
    DIRS[f'a{i+1}, b{i+1}'] = e((3 + i, 1))
    DIRS[f'c{i+1}'] = e((i, 1), (3 + i, -1))
    DIRS[f'alpha{i+1}, beta{i+1}'] = e((6 + i, 1))
    DIRS[f'gamma{i+1}'] = e((i, -1), (6 + i, -1))

def certify(cert):
    mp.dps=90;iv.dps=90;Jet.ctx=iv
    image=iv.mpf(cert['image_radius_upper'].strip('[]').split(',')[1])
    if not cert['strict_inclusion'] or not image.b<iv.mpf('3.1e-30').a:raise ArithmeticError('Fresh root image does not fit the rank box.')
    J=jac('3.1e-30',cert);out={}
    for name,b in DIRS.items():
        det,piv=kalman(J,b)
        good=det is not None and (det.a>0 or det.b<0)
        if not good:raise ArithmeticError('Unresolved Kalman rank: '+name)
        out[name]=dict(direction=b,rank9=bool(good),scaled_determinant=str(det),pivots=len(piv))
    det,_=kalman(jac('1e-20',cert),DIRS['alpha1, beta1'])
    if not (det.a>0 or det.b<0):raise ArithmeticError('Full-box alpha1 rank unresolved.')
    return dict(inputs=out,all_eighteen_rates=True,alpha1_full_box=str(det),scope='Fresh directed intervals over the freshly certified generalized-Hopf root box; rank is not a conditioning or finite-witness switching certificate.')