KiHoLee
17d23fa76a
check_family_enum.py measures the attack the manuscript now states in Section III-A: the winning correlation is an index-free verifier, so ranking the 63 non-constant Walsh rows by mean winning correlation recovers the user set from one frame in 0.905 of 200 trials at 10 dB and from four frames in 0.990, using nothing outside the stated threat model. Under the invariance refresh it recovers it in none, because the entry permutation relabels the codebook the adversary must align against. V8 and V9 read the trained codebook through main_model(), which retrains on every call, and a codebook trained on CUDA is not the one trained on CPU. The shipped verify_math.csv therefore read PASS here and FAIL for anyone running this package without a GPU. model_main.pt is 7 KB and fixes the codebook, which is what both checks are about; delete it to retrain. V1-V11 now pass on both. New checks: V10, the format-matched OMA reference Section VI-B quotes, and V11, the closed-form against Monte Carlo comparison the manuscript claimed and never stored. V3a's bias-linearity result was computed and printed but never written to the CSV, so the one linearity claim the paper quotes was the one this package could not show. check_consistency.py gains 21 assertions, covering five data files that no assertion read (users, csi, semantic, cov_attack, sec_jam) and the trend claims it structurally could not see, since it compared values and not shapes. README: the figure map named stages that do not write the artifacts they list, so following it did not reproduce Figs. 4 and 6; the reproduction block was five scripts short; and the refresh numbers were from a superseded run (nearly three, 15.0 to 64.8 bits) against the manuscript's 2.3 and 23.8 to 364.6.
406 lines
16 KiB
Python
406 lines
16 KiB
Python
"""Numerical verification of the closed forms in Sections IV and V of
|
|
paper 11 (mask-as-key encryption and jamming robustness).
|
|
|
|
Every claim that enters the manuscript is checked here against Monte
|
|
Carlo, with a PASS/FAIL verdict and the achieved agreement level printed.
|
|
Real-vector convention, dimension d, U users, codebook of V unit-norm
|
|
codewords, per-user masks with zero-mean entries normalized to
|
|
||m||^2 = d (so E[m_k^2] = 1).
|
|
|
|
Notation matches the tex:
|
|
y = (1/c) sum_u e_{s_u} .* m_u transmit frame
|
|
legit score z_{u,i} = r_u^T (e_i .* m_u), r_u = y + n_u/h_u
|
|
eve score zE_{u,i} = (y_E/h_E)^T (e_i .* mtil_u)
|
|
jammer adds h_J sqrt(rho) w to the victim observation
|
|
|
|
Run on CPU (NumPy); no training involved, pure algebra checks.
|
|
"""
|
|
from __future__ import annotations
|
|
import numpy as np
|
|
|
|
from pathlib import Path
|
|
|
|
RNG = np.random.default_rng(2026)
|
|
D, U, V = 64, 4, 256
|
|
|
|
CKPT = Path(__file__).resolve().parent.parent / "data" / "model_main.pt"
|
|
|
|
|
|
def cached_main_model():
|
|
"""The trained main-configuration model, from a checkpoint.
|
|
|
|
V8 and V9 read the trained codebook. Retraining it reproduces only
|
|
on the device that trained it, so a CPU run of the released package
|
|
disagreed with the shipped numbers. The checkpoint fixes the
|
|
codebook, which is what both checks are about; delete it to retrain.
|
|
"""
|
|
import torch
|
|
from exp_full import main_model
|
|
m = main_model()
|
|
if CKPT.exists():
|
|
m.load_state_dict(torch.load(CKPT, map_location="cpu"))
|
|
else:
|
|
torch.save({k: v.cpu() for k, v in m.state_dict().items()}, CKPT)
|
|
return m
|
|
|
|
|
|
|
|
def unit_codebook(V, d, rng):
|
|
E = rng.standard_normal((V, d))
|
|
return E / np.linalg.norm(E, axis=1, keepdims=True)
|
|
|
|
|
|
def masks(U, d, rng):
|
|
"""Zero-mean entries normalized so ||m_u||^2 = d."""
|
|
M = rng.standard_normal((U, d))
|
|
return M / np.linalg.norm(M, axis=1, keepdims=True) * np.sqrt(d)
|
|
|
|
|
|
ROWS = [] # (tag, claim, emp, err, tol, verdict)
|
|
|
|
|
|
def report(tag, claim, emp, tol, extra=""):
|
|
err = abs(claim - emp)
|
|
ok = err <= tol
|
|
ROWS.append((tag, claim, emp, err, tol, "PASS" if ok else "FAIL"))
|
|
print(f"[{'PASS' if ok else 'FAIL'}] {tag}: claim={claim:.5g} "
|
|
f"emp={emp:.5g} |err|={err:.2g} tol={tol:g} {extra}")
|
|
return ok
|
|
|
|
|
|
def v1_legit_self_alignment():
|
|
"""Claim: E[e_s^T diag(m^2) e_s] = 1 (signal self-correlation)."""
|
|
E = unit_codebook(V, D, RNG)
|
|
vals = []
|
|
for _ in range(4000):
|
|
m = masks(1, D, RNG)[0]
|
|
s = RNG.integers(V)
|
|
vals.append(float((E[s] ** 2) @ (m ** 2)))
|
|
return report("V1 legit self-alignment", 1.0, float(np.mean(vals)), 2e-2)
|
|
|
|
|
|
def v2_eve_uninformed():
|
|
"""Claim: with an independent substitute mask, the eavesdropper's
|
|
correct-index correlation has the same mean as any wrong index, so
|
|
the mean advantage is zero and the eavesdropper SER = (V-1)/V,
|
|
independent of SNR."""
|
|
E = unit_codebook(V, D, RNG)
|
|
# mean advantage of the true index over the wrong indices, noiseless
|
|
adv = []
|
|
ser_by_snr = {}
|
|
for snr_db in [0.0, 10.0, 20.0, 80.0]: # 80 dB stands in for noiseless
|
|
sigma = np.sqrt(1.0 / (D * 10 ** (snr_db / 10.0)))
|
|
err = 0
|
|
trials = 6000
|
|
for _ in range(trials):
|
|
s = RNG.integers(V, size=U)
|
|
M = masks(U, D, RNG)
|
|
c = 1.0 # scale-invariant for argmax
|
|
y = np.zeros(D)
|
|
for u in range(U):
|
|
y += E[s[u]] * M[u]
|
|
y /= np.sqrt(U) # any fixed scale
|
|
# eavesdropper targets user 0 with an independent wrong mask
|
|
mtil = masks(1, D, RNG)[0]
|
|
hE = np.sqrt(-np.log(RNG.random()))
|
|
rE = y + (sigma / hE) * RNG.standard_normal(D)
|
|
scores = (E * mtil) @ rE # (V,)
|
|
if snr_db == 80.0:
|
|
adv.append(scores[s[0]] - scores.mean())
|
|
if scores.argmax() != s[0]:
|
|
err += 1
|
|
ser_by_snr[snr_db] = err / trials
|
|
chance = (V - 1) / V
|
|
ok1 = report("V2a eve mean advantage", 0.0, float(np.mean(adv)), 3e-3)
|
|
ok2 = True
|
|
for snr_db, ser in ser_by_snr.items():
|
|
tag = f"V2b eve SER @ {int(snr_db)}dB"
|
|
ok2 &= report(tag, chance, ser, 1.5e-2)
|
|
return ok1 and ok2
|
|
|
|
|
|
def v3_leakage_vs_correlation():
|
|
"""Claim: if the substitute mask has normalized correlation
|
|
rho = <m,mtil>/d with the true mask, the eavesdropper's true-index
|
|
bias grows linearly in rho; independent random masks give
|
|
E|rho| = O(1/sqrt(d)); orthogonal masks give rho = 0."""
|
|
E = unit_codebook(V, D, RNG)
|
|
# (a) bias vs prescribed rho
|
|
slopes = []
|
|
for rho in [0.0, 0.25, 0.5, 0.75, 1.0]:
|
|
bias = []
|
|
for _ in range(3000):
|
|
m = masks(1, D, RNG)[0]
|
|
mp = masks(1, D, RNG)[0]
|
|
mp = mp - (mp @ m) / (m @ m) * m # orthogonalize
|
|
mp = mp / np.linalg.norm(mp) * np.sqrt(D)
|
|
mtil = rho * m + np.sqrt(1 - rho ** 2) * mp
|
|
s = RNG.integers(V)
|
|
# noiseless single-user useful alignment for the true index
|
|
bias.append(float((E[s] ** 2) @ (m * mtil)))
|
|
slopes.append((rho, float(np.mean(bias))))
|
|
# claim: bias(rho) = rho * bias(1); check linearity
|
|
b1 = slopes[-1][1]
|
|
lin_ok = all(abs(b - rho * b1) <= 3e-2 for rho, b in slopes)
|
|
print(f"[{'PASS' if lin_ok else 'FAIL'}] V3a bias linear in rho: "
|
|
+ ", ".join(f"rho={r:.2f}->{b:.3f}" for r, b in slopes))
|
|
ROWS.append(("V3a bias slope in kappa", "1.0", "%.4f" % b1,
|
|
"%.4f" % abs(1.0 - b1), "0.03",
|
|
"PASS" if lin_ok else "FAIL"))
|
|
# (b) random independent mask correlation: E|corr| = sqrt(2/(pi d))
|
|
# (the folded-normal mean of a N(0, 1/d) variable)
|
|
corrs = []
|
|
for _ in range(5000):
|
|
m = masks(1, D, RNG)[0]
|
|
mt = masks(1, D, RNG)[0]
|
|
corrs.append(abs((m @ mt) / D))
|
|
emp = float(np.mean(corrs))
|
|
claim = float(np.sqrt(2.0 / (np.pi * D)))
|
|
ok_b = report("V3b random mask E|corr|", claim, emp, 0.03 * claim)
|
|
return lin_ok and ok_b
|
|
|
|
|
|
def v4_blind_jammer_spread():
|
|
"""Claim: a mask-blind jammer (w independent of m_u) contributes a
|
|
zero-mean term to every candidate score with variance
|
|
(hJ^2/hu^2) rho * sum_k w_k^2 e_{i,k}^2, i.e. it is spread with no
|
|
systematic bias toward any index."""
|
|
E = unit_codebook(V, D, RNG)
|
|
rho = 1.0
|
|
# (a) structural claim: the mask projection g_i = (w .* m)^T e_i is
|
|
# zero-mean, decoupled from the positive gain factor hJ/hu.
|
|
g, var_emp, var_cl = [], [], []
|
|
for _ in range(20000):
|
|
m = masks(1, D, RNG)[0]
|
|
w = RNG.standard_normal(D); w /= np.linalg.norm(w)
|
|
i = RNG.integers(V)
|
|
gi = (w * m) @ E[i]
|
|
g.append(gi)
|
|
var_emp.append(gi ** 2)
|
|
var_cl.append(np.sum(w ** 2 * E[i] ** 2))
|
|
ok1 = report("V4a blind jammer projection mean", 0.0,
|
|
float(np.mean(g)), 3e-3)
|
|
# (b) variance of the projection matches sum_k w_k^2 e_{i,k}^2; the
|
|
# full contribution scales this by (hJ^2/hu^2) rho.
|
|
ok2 = report("V4b blind jammer projection variance",
|
|
float(np.mean(var_cl)), float(np.mean(var_emp)),
|
|
0.01 * float(np.mean(var_cl)))
|
|
return ok1 and ok2
|
|
|
|
|
|
def v5_matched_jammer_concentrates():
|
|
"""Claim: a mask-matched jammer aligned with the victim key for a
|
|
target index t creates a bias of order one on index t, while the
|
|
blind-jammer projection has zero mean and RMS of order 1/sqrt(d).
|
|
The physically meaningful separation is matched bias over blind RMS,
|
|
which is sqrt(d) (the sample mean of the blind bias estimates zero
|
|
and is pure Monte Carlo noise, so it is NOT a valid denominator)."""
|
|
E = unit_codebook(V, D, RNG)
|
|
bias_matched, blind_sq = [], []
|
|
for _ in range(3000):
|
|
m = masks(1, D, RNG)[0]
|
|
t = RNG.integers(V)
|
|
wm = E[t] * m; wm /= np.linalg.norm(wm) # matched (needs m)
|
|
wb = RNG.standard_normal(D); wb /= np.linalg.norm(wb) # blind
|
|
bias_matched.append(float((wm * m) @ E[t]))
|
|
blind_sq.append(float(((wb * m) @ E[t]) ** 2))
|
|
bm = float(np.mean(bias_matched))
|
|
brms = float(np.sqrt(np.mean(blind_sq)))
|
|
ratio = bm / brms
|
|
ok = abs(ratio - np.sqrt(D)) <= 0.25 * np.sqrt(D) and bm > 0.9
|
|
print(f"[{'PASS' if ok else 'FAIL'}] V5 matched bias / blind RMS: "
|
|
f"matched={bm:.3f} blind_rms={brms:.4f} ratio={ratio:.1f} "
|
|
f"(claim sqrt(d)={np.sqrt(D):.1f})")
|
|
ROWS.append(("V5 matched bias over blind RMS",
|
|
"%.1f" % np.sqrt(D), "%.1f" % ratio,
|
|
"%.2f" % abs(ratio - np.sqrt(D)), "1.0",
|
|
"PASS" if ok else "FAIL"))
|
|
return ok
|
|
|
|
|
|
def v6_coded_oma_outage():
|
|
"""The genie reference the manuscript concedes: an OMA user coding at
|
|
the Rayleigh outage limit of its own allocation.
|
|
|
|
The user owns L = d/U real dimensions, so d/(2U) complex uses, and
|
|
carries log2(V) bits. Its per-dimension SNR is the frame SNR, the
|
|
same convention snr_to_sigma2 sets. Outage is the probability that
|
|
the instantaneous mutual information falls below that rate."""
|
|
import math
|
|
U, d, V = 4, 256, 65536
|
|
for snr_db in (10.0,):
|
|
g = 10.0 ** (snr_db / 10.0)
|
|
uses = d / (2.0 * U) # complex channel uses
|
|
rate = math.log2(V) / uses # bits per complex use
|
|
thr = (2.0 ** rate - 1.0) / g # |h|^2 threshold
|
|
pout = 1.0 - math.exp(-thr) # |h|^2 ~ Exp(1)
|
|
print(f"V6 coded-OMA outage @ {snr_db:.0f} dB: {pout:.4f} "
|
|
f"(rate {rate:.3f} bit/use)")
|
|
ROWS.append(("V6 coded-OMA outage @ %.0f dB" % snr_db,
|
|
"%.4f" % pout, "%.6f" % pout, "0", "0", "REFERENCE"))
|
|
return True
|
|
|
|
|
|
def v7_symbolic_identities():
|
|
"""Symbolic verification of the three algebraic identities. Monte
|
|
Carlo cannot check an identity, only an instance of it."""
|
|
try:
|
|
import sympy as sp
|
|
except ImportError:
|
|
print("V7 symbolic: sympy not installed, SKIPPED")
|
|
ROWS.append(("V7 symbolic identities", "-", "-", "-", "-", "SKIPPED"))
|
|
return True
|
|
ok = True
|
|
|
|
# refresh entropy: L sign bits, an entry permutation, a user permutation
|
|
L, Uu = 64, 4
|
|
ent = sp.Integer(L) + sp.log(sp.factorial(L), 2) + sp.log(sp.factorial(Uu), 2)
|
|
ok &= abs(float(ent) - 364.6) < 0.05
|
|
print("V7a refresh entropy %.3f bits/block" % float(ent))
|
|
|
|
# fixed-key entropy: ordered choices of U of the L-1 non-constant rows
|
|
fixed = sp.log(sp.factorial(L - 1) / sp.factorial(L - 1 - Uu), 2)
|
|
ok &= abs(float(fixed) - 23.8) < 0.05
|
|
ok &= sp.factorial(L - 1) / sp.factorial(L - 1 - Uu) == 14295960
|
|
print("V7b fixed-key entropy %.3f bits, %d choices"
|
|
% (float(fixed), sp.factorial(L - 1) / sp.factorial(L - 1 - Uu)))
|
|
|
|
# the termwise Hadamard identity behind eq:hadamard
|
|
n = 4
|
|
e = sp.Matrix(sp.symbols("e1:%d" % (n + 1)))
|
|
m = sp.Matrix(sp.symbols("m1:%d" % (n + 1)))
|
|
f = sp.Matrix(sp.symbols("f1:%d" % (n + 1)))
|
|
lhs = sum((sp.matrix_multiply_elementwise(e, m))[k] *
|
|
(sp.matrix_multiply_elementwise(f, m))[k] for k in range(n))
|
|
rhs = (e.T * sp.diag(*[m[k] ** 2 for k in range(n)]) * f)[0]
|
|
ok &= sp.simplify(lhs - rhs) == 0
|
|
print("V7c (e*m).(f*m) == e^T diag(m^2) f :",
|
|
sp.simplify(lhs - rhs) == 0)
|
|
|
|
ROWS.append(("V7 symbolic identities", "exact", "exact", "0", "0",
|
|
"PASS" if ok else "FAIL"))
|
|
return ok
|
|
|
|
|
|
def v8_cross_period_terms():
|
|
"""Proposition 2 keeps only the diagonal of the jammer projection.
|
|
The periodic key makes entries one period apart identical, so the
|
|
cross-period terms do not vanish termwise. They are zero mean over
|
|
the codebook, which is the claim the proof rests on."""
|
|
import math
|
|
import torch
|
|
m = cached_main_model()
|
|
Bn = m.unit_codebook().detach().cpu()
|
|
pat = m.masks().detach().cpu()[0]
|
|
L, P, d = m.L, m.P, m.d
|
|
# a generator of its own, seeded after the model is built: seeding the
|
|
# global one first leaves the draw dependent on how main_model consumed
|
|
# it, which moved this number between runs
|
|
g = torch.Generator().manual_seed(7)
|
|
rel = []
|
|
for _ in range(20000):
|
|
w = torch.randn(d, generator=g)
|
|
w /= w.norm()
|
|
i = torch.randint(m.vu, (P,), generator=g)
|
|
e = (Bn[i] / math.sqrt(P)).reshape(-1)
|
|
a = (w * e).reshape(P, L) * pat[None, :]
|
|
diag = float((a ** 2).sum())
|
|
rel.append((float((a.sum(0) ** 2).sum()) - diag) / diag)
|
|
mean = sum(rel) / len(rel)
|
|
ok = abs(mean) < 0.0005
|
|
print("V8 cross-period remainder, mean %+.4f of the retained term" % mean)
|
|
ROWS.append(("V8 cross-period remainder", "0.0", "%.6f" % mean,
|
|
"%.6f" % abs(mean), "0.0005", "PASS" if ok else "FAIL"))
|
|
return ok
|
|
|
|
|
|
def v9_score_variance_ratio():
|
|
"""Per-digit score-variance ratio of Proposition 1's proof: the mean
|
|
of sum_j e_j^4 / sum_j e_j^2 e'_j^2 over ordered codeword pairs of
|
|
the trained unit codebook, quoted as 2.8 in the manuscript."""
|
|
import torch
|
|
m = cached_main_model()
|
|
B = m.unit_codebook().detach().cpu().double()
|
|
B = B / B.norm(dim=1, keepdim=True)
|
|
n = B.shape[0]
|
|
num = (B ** 4).sum(1)
|
|
ratios = []
|
|
for i in range(n):
|
|
for j in range(n):
|
|
if i != j:
|
|
ratios.append(float(num[i] / (B[i] ** 2 * B[j] ** 2).sum()))
|
|
mean = sum(ratios) / len(ratios)
|
|
ok = abs(mean - 2.8) < 0.05
|
|
print("V9 per-digit score-variance ratio: %.3f (quoted 2.8)" % mean)
|
|
ROWS.append(("V9 score-variance ratio", "2.8", "%.4f" % mean,
|
|
"%.4f" % abs(mean - 2.8), "0.05", "PASS" if ok else "FAIL"))
|
|
return ok
|
|
|
|
|
|
|
|
def v10_format_matched_oma():
|
|
"""The format-matched OMA reference of Section VI-B.
|
|
|
|
The binary reference spends 16 of its 64 exclusive dimensions on
|
|
antipodal bits. The same allocation spent the way the proposed
|
|
scheme spends it, P=4 sixteen-ary orthogonal decisions over 16
|
|
dimensions each, is the comparison a reviewer will ask for."""
|
|
from sse_lib import oma_ser_orth
|
|
from exp_full import oma_ser_keylen
|
|
val = oma_ser_orth([10.0])[0]
|
|
binary = oma_ser_keylen(64, 10.0)
|
|
ok = abs(val - 0.055) < 0.001
|
|
print("V10 format-matched OMA at 10 dB: %.5f (binary %.5f)"
|
|
% (val, binary))
|
|
ROWS.append(("V10 format-matched OMA at 10 dB", "0.055", "%.5f" % val,
|
|
"%.5f" % abs(val - 0.055), "0.001", "PASS" if ok else "FAIL"))
|
|
return ok
|
|
|
|
|
|
def v11_oma_closed_form_vs_mc():
|
|
"""The manuscript says the OMA closed form agrees with Monte Carlo
|
|
to within one percent. That check had no stored artifact."""
|
|
from sse_lib import oma_ser, oma_ser_mc
|
|
cf = oma_ser([16.0])[0]
|
|
mc = oma_ser_mc([16.0], frames=2_000_000)[0]
|
|
rel = abs(cf - mc) / mc
|
|
ok = rel < 0.01
|
|
print("V11 OMA closed form %.6f vs Monte Carlo %.6f (%.2f%%)"
|
|
% (cf, mc, 100 * rel))
|
|
ROWS.append(("V11 OMA closed form vs Monte Carlo", "%.6f" % mc,
|
|
"%.6f" % cf, "%.4f" % rel, "0.01", "PASS" if ok else "FAIL"))
|
|
return ok
|
|
|
|
|
|
def main():
|
|
print(f"config d={D} U={U} V={V}\n")
|
|
results = {
|
|
"V1": v1_legit_self_alignment(),
|
|
"V2": v2_eve_uninformed(),
|
|
"V3": v3_leakage_vs_correlation(),
|
|
"V4": v4_blind_jammer_spread(),
|
|
"V5": v5_matched_jammer_concentrates(),
|
|
"V6": v6_coded_oma_outage(),
|
|
"V7": v7_symbolic_identities(),
|
|
"V8": v8_cross_period_terms(),
|
|
"V9": v9_score_variance_ratio(),
|
|
"V10": v10_format_matched_oma(),
|
|
"V11": v11_oma_closed_form_vs_mc(),
|
|
}
|
|
print("\nsummary:", {k: ("PASS" if v else "FAIL") for k, v in results.items()})
|
|
print("ALL PASS" if all(results.values()) else "SOME FAILED")
|
|
# stored artifact so every quoted verification number has a raw file
|
|
import csv as _csv
|
|
from pathlib import Path as _Path
|
|
data = _Path(__file__).resolve().parents[1] / "data"
|
|
with open(data / "verify_math.csv", "w", newline="") as f:
|
|
w = _csv.writer(f)
|
|
w.writerow(["check", "claim", "empirical", "abs_err", "tol",
|
|
"verdict"])
|
|
w.writerows(ROWS)
|
|
print("[csv]", data / "verify_math.csv")
|
|
|
|
|
|
if __name__ == "__main__":
|
|
main()
|