#!/usr/bin/env python3
"""Comprobaciones del experimento sha_exp.py, hechas por el verificador el 3 de noviembre de 2026.

Uso (con numpy y scikit-learn, desde un directorio de trabajo fuera del repositorio):
    python -I sha_exp_verif.py vectores     # NIST, hash real y recorte de rondas contra otra implementación
    python -I sha_exp_verif.py solape       # nonces repetidos entre entrenamiento y prueba
    python -I sha_exp_verif.py base         # línea base de clase mayoritaria y señal por registro
    python -I sha_exp_verif.py diferencias  # test clásico sin IA: cambiar un bit del nonce
    python -I sha_exp_verif.py dedicado     # red y árboles solo con los registros menos mezclados (d y h)

Los registros se nombran a..h como en FIPS 180-4. Tras r rondas, d = a(r-3) y h = e(r-3) más una
constante, de modo que un bit de la ronda r-3 sigue visible en la salida de la ronda r.
"""
import hashlib, importlib.util, json, os, struct, sys, time
import numpy as np

sys.dont_write_bytecode = True  # -I ignora PYTHONDONTWRITEBYTECODE; así no se deja __pycache__ en el repositorio
_p = os.path.join(os.path.dirname(os.path.abspath(__file__)), 'sha_exp.py')
_s = importlib.util.spec_from_file_location('sha_exp', _p)
sx = importlib.util.module_from_spec(_s)
_s.loader.exec_module(sx)

K = [int(k) for k in sx.K]
IV = [int(k) for k in sx.IV]
M = 2**32 - 1
REG = 'abcdefgh'


def rr(x, n):
    return ((x >> n) | (x << (32 - n))) & M


def comp(st, blk, r=64):
    """Compresión de SHA-256 en Python puro, sin numpy, con r rondas."""
    W = list(struct.unpack('>16I', blk))
    for t in range(16, 64):
        W.append((W[t-16] + (rr(W[t-15], 7) ^ rr(W[t-15], 18) ^ (W[t-15] >> 3)) + W[t-7]
                  + (rr(W[t-2], 17) ^ rr(W[t-2], 19) ^ (W[t-2] >> 10))) & M)
    a, b, c, d, e, f, g, h = st
    for t in range(r):
        t1 = (h + (rr(e, 6) ^ rr(e, 11) ^ rr(e, 25)) + ((e & f) ^ ((~e) & M & g)) + K[t] + W[t]) & M
        t2 = ((rr(a, 2) ^ rr(a, 13) ^ rr(a, 22)) + ((a & b) ^ (a & c) ^ (b & c))) & M
        h, g, f, e, d, c, b, a = g, f, e, (d + t1) & M, c, b, a, (t1 + t2) & M
    return [(x + y) & M for x, y in zip(st, [a, b, c, d, e, f, g, h])]


def pad(m):
    n = len(m)
    m = m + b'\x80'
    m += b'\0' * ((56 - len(m)) % 64)
    return m + struct.pack('>Q', n * 8)


def vectores():
    casos = ((b'abc', 'ba7816bf8f01cfea414140de5dae2223b00361a396177a9cb410ff61f20015ad'),
             (b'', 'e3b0c44298fc1c149afbf4c8996fb92427ae41e4649b934ca495991b7852b855'),
             (b'abcdbcdecdefdefgefghfghighijhijkijkljklmklmnlmnomnopnopq',
              '248d6a61d20638b8e5c026930c3e6039a33ce45964ff2167f6ecedd419db06c1'))
    for msg, esperado in casos:
        p = pad(msg)
        st = IV
        st2 = [np.uint32(x) for x in IV]
        for i in range(0, len(p), 64):
            st = comp(st, p[i:i+64])
            o = sx.compress(st2, [np.uint32(x) for x in struct.unpack('>16I', p[i:i+64])])
            st2 = [np.uint32(x[0]) for x in o]
        d1 = b''.join(x.to_bytes(4, 'big') for x in st).hex()
        d2 = b''.join(int(x).to_bytes(4, 'big') for x in st2).hex()
        assert d1 == d2 == esperado == hashlib.sha256(msg).hexdigest()
    print('vectores FIPS 180-4: correctos con las dos implementaciones')
    mid = comp(IV, sx.HDR[:64])
    assert [int(x) for x in sx.MID] == mid
    ns = np.random.default_rng(1).integers(0, 2**32, 200, dtype=np.uint32)
    for r in range(1, 65):
        o = sx.block2(ns, rounds=r)
        for i, n in enumerate(ns):
            blk = sx.HDR[64:76] + int(n).to_bytes(4, 'little') + b'\x80' + b'\0' * 39 + struct.pack('>Q', 640)
            assert [int(x[i]) for x in o] == comp(mid, blk, r), (r, i)
    print('recorte de rondas 1 a 64: coincide con la implementación en Python puro en 200 nonces')
    for r in range(1, 6):
        o = sx.block2(ns, rounds=r)
        print(' ronda', r, 'registros que dependen del nonce:', [REG[i] for i, x in enumerate(o) if np.unique(x).size > 1])


def solape():
    g = np.random.default_rng(974323)
    for ntr, nte, nom in ((300000, 50000, 'A'), (2000000, 1000000, 'B')):
        x = g.integers(0, 2**32, ntr + nte, dtype=np.uint32)
        c = int(np.isin(x[ntr:], x[:ntr]).sum())
        print(f'experimento {nom}: {c} de {nte} nonces de prueba están también en entrenamiento ({c/nte:.5%})')


def base():
    rng = np.random.default_rng(974323)
    ntr, nte = 300000, 50000
    um = 3.9 * .5 / np.sqrt(nte)
    for r in (4, 5, 6, 7, 8, 9, 10):
        nonce = rng.integers(0, 2**32, ntr + nte, dtype=np.uint32)
        o = sx.block2(nonce, rounds=r)
        Y = np.concatenate([sx.bits(x) for x in o], axis=1)
        X = sx.bits(nonce) * 2 - 1
        frac = Y[:ntr].mean(0)
        var = (frac > .02) & (frac < .98)
        idx = np.where(var)[0]
        m = sx.Red(32, int(var.sum())).fit(X[:ntr], Y[:ntr][:, var])
        acc = ((m.adelante(X[ntr:])[2] > 0) == (Y[ntr:][:, var] > .5)).mean(0)
        maj = np.maximum(frac[var], 1 - frac[var])
        sen = acc > .5 + um
        print(dict(r=r, variables=int(var.sum()), sobre_50=int(sen.sum()), sobre_mayoritaria=int((acc > maj + um).sum()),
                   mayoritaria_max=round(float(maj.max()), 4), por_registro={REG[i]: int(sen[(idx // 32) == i].sum()) for i in range(8)}), flush=True)


def diferencias(n=100000):
    base_n = np.random.default_rng(7).integers(0, 2**32, n, dtype=np.uint32)
    res = []
    for r in list(range(4, 25)) + [28, 32, 48, 64]:
        o = sx.block2(base_n, rounds=r)
        dmax, celdas = 0, 0
        for b in range(32):
            o2 = sx.block2(base_n ^ (np.uint32(1) << np.uint32(b)), rounds=r)
            p = np.concatenate([sx.bits(x ^ y) for x, y in zip(o, o2)], axis=1).mean(0)
            d = np.abs(p - .5)
            dmax = max(dmax, float(d.max()))
            celdas += int((d > 5 * .5 / np.sqrt(n)).sum())
        rec = dict(rondas=r, desvio_max=round(dmax, 4), celdas_con_senal=celdas, de=32 * 256)
        res.append(rec)
        print(rec, flush=True)
    json.dump(res, open('diferencias.json', 'w'), indent=1)


def dedicado():
    from sklearn.ensemble import HistGradientBoostingClassifier as G
    rng = np.random.default_rng(20261103)
    ntr, nte = 1500000, 200000
    um = 3.9 * .5 / np.sqrt(nte)
    for r in (9, 10):
        nonce = rng.integers(0, 2**32, ntr + nte, dtype=np.uint32)
        o = sx.block2(nonce, rounds=r)
        Y = np.concatenate([sx.bits(o[3]), sx.bits(o[7])], axis=1)
        X = sx.bits(nonce) * 2 - 1
        m = sx.Red(32, 64, sem=3).fit(X[:ntr], Y[:ntr], ep=12, lote=512)
        acc = ((m.adelante(X[ntr:])[2] > 0) == (Y[ntr:] > .5)).mean(0)
        print(dict(r=r, modelo='red, 64 salidas, 1,5 M', con_senal=int((acc > .5 + um).sum()), acierto_max=round(float(acc.max()), 4)), flush=True)
        acct = []
        modelos = []
        for j in range(64):
            g = G(max_iter=150, learning_rate=.1, max_leaf_nodes=63, random_state=1).fit(X[:600000] > 0, Y[:600000, j])
            acct.append(float((g.predict(X[ntr:] > 0) == Y[ntr:, j]).mean()))
            modelos.append(g)
        acct = np.array(acct)
        j = int(acct.argmax())
        nom = ('d' if j < 32 else 'h') + str(j % 32)
        print(dict(r=r, modelo='árboles por bit, 600 k', con_senal=int((acct > .5 + um).sum()), mejor_bit=nom, acierto_max=round(float(acct.max()), 4)), flush=True)
        # Se confirma el mejor bit con nonces nuevos, que no han servido para elegirlo
        nv = np.random.default_rng(r * 1000 + 1).integers(0, 2**32, 2000000, dtype=np.uint32)
        o2 = sx.block2(nv, rounds=r)
        y2 = sx.bits(o2[3] if j < 32 else o2[7])[:, j % 32]
        a2 = float((modelos[j].predict(sx.bits(nv) > .5) == y2).mean())
        print(dict(r=r, confirmacion=nom, nonces_nuevos=len(nv), acierto=round(a2, 5), desviaciones=round((a2 - .5) / (.5 / np.sqrt(len(nv))), 1)), flush=True)


if __name__ == '__main__':
    {'vectores': vectores, 'solape': solape, 'base': base, 'diferencias': diferencias, 'dedicado': dedicado}[sys.argv[1]]()
