diff --git a/code/diag_interference.py b/code/diag_interference.py new file mode 100644 index 0000000..d37cc64 --- /dev/null +++ b/code/diag_interference.py @@ -0,0 +1,88 @@ +# -*- coding: utf-8 -*- +"""Where does the legitimate advantage over OMA go? + +An ideal M-ary receiver at the main configuration should reach 0.199 at +10 dB against the 0.275 of resource-matched OMA, a factor of 1.38, while +the system measures 0.257, a factor of 1.07. This script splits the +shortfall into its two causes: residual multi-user interference, which +orthogonal keys do not remove because masking is elementwise, and the +distance the trained unit codebook falls short of an orthogonal set. +""" +import math +import sys +from pathlib import Path + +import torch + +sys.path.insert(0, str(Path(__file__).resolve().parent)) +import sse_lib as L +from sse_lib import DEVICE, snr_to_sigma2, rayleigh_gain +from exp_full import main_model + +SNR_DB = 10.0 +FRAMES = 400_000 +CH = 40_000 + + +def ser(m, snr_db, frames, solo=False): + """SER of user 0. With solo=True the other users transmit nothing, + while the power normalizer c is left at its four-user value so that + user 0 keeps exactly the energy it has in the real system.""" + tot = wrong = 0 + with torch.no_grad(): + while tot < frames: + n = min(CH, frames - tot) + dig = torch.randint(m.vu, (n, m.users, m.P), device=DEVICE) + Bn = m.unit_codebook() + e = Bn[dig] / math.sqrt(m.P) + mk = m.masks() + x = e * mk[None, :, None, :] + if solo: + x = x[:, :1] + y = x.sum(dim=1) / m.c + h = rayleigh_gain((n, 1), device=DEVICE) + sig = snr_to_sigma2(torch.full((n,), snr_db), m.d).to(DEVICE).sqrt() + rx = h[:, :, None, None] * y[:, None] \ + + sig[:, None, None, None] * torch.randn(n, 1, m.P, m.L, + device=DEVICE) + r = rx / h[:, :, None, None].clamp_min(1e-6) + cand = Bn[None, :, :] * mk[:1, None, :] + sc = torch.einsum("nupl,uvl->nupv", r, cand) + bad = (sc.argmax(-1)[:, 0] != dig[:, 0]).any(dim=-1) + wrong += int(bad.sum()) + tot += n + return wrong / tot + + +def main(): + m = main_model() + with torch.no_grad(): + Bn = m.unit_codebook() + G = Bn @ Bn.T + off = G - torch.diag(torch.diag(G)) + mk = m.masks() + Gm = mk @ mk.T / m.L + offm = Gm - torch.diag(torch.diag(Gm)) + + print("main configuration: d=%d P=%d L=%d Vu=%d U=%d" + % (m.d, m.P, m.L, m.vu, m.users)) + print("key cross-correlation, max |off-diagonal| : %.2e" + % offm.abs().max()) + print("codebook Gram, max |off-diagonal| : %.4f" + % off.abs().max()) + print("codebook Gram, rms off-diagonal : %.4f" + % off.pow(2).sum().div(m.vu * (m.vu - 1)).sqrt()) + print("(an orthogonal set of %d codewords in %d dims would read 0)" + % (m.vu, m.L)) + print() + four = ser(m, SNR_DB, FRAMES, solo=False) + solo = ser(m, SNR_DB, FRAMES, solo=True) + print("user-0 SER, all four users transmitting : %.4f" % four) + print("user-0 SER, other users silent : %.4f" % solo) + print("OMA, resource matched (closed form) : %.4f" + % L.oma_ser([SNR_DB])[0]) + print("ideal 16-ary orthogonal (separate MC) : 0.1986") + + +if __name__ == "__main__": + main() diff --git a/code/diag_orthobook.py b/code/diag_orthobook.py new file mode 100644 index 0000000..adf4b7d --- /dev/null +++ b/code/diag_orthobook.py @@ -0,0 +1,105 @@ +# -*- coding: utf-8 -*- +"""Does an orthogonal unit codebook recover the shortfall? + +diag_interference shows the gap to the ideal M-ary receiver is not +multi-user interference but the geometry of the trained unit codebook, +whose Gram matrix carries a root-mean-square off-diagonal of 0.45 where +an orthogonal set would carry zero. Vu = L = 16 admits an exactly +orthogonal set, so this measures what installing one buys. + +Two orthogonal sets are tried, because the choice is not free. The +Walsh-Hadamard set collides with the keys: the rows are closed under the +elementwise product, so masking a Hadamard codeword by a Hadamard key +returns another Hadamard codeword and every user ends up with the same +candidate set. A random orthogonal set carries no such group structure, +and masking by a unit-modulus key preserves its orthogonality exactly. +""" +import sys +from pathlib import Path + +import torch + +sys.path.insert(0, str(Path(__file__).resolve().parent)) +import sse_lib as L +from sse_lib import DEVICE, SSE +from exp_full import hadamard, base_keys +from diag_interference import ser + +SNR = [0.0, 10.0, 20.0] +FRAMES = 400_000 + + +def fixed_model(B, P=4, vu=16, d=64, U=4): + L.set_seed(1) + m = SSE(P=P, vu=vu, d=d, users=U).to(DEVICE) + with torch.no_grad(): + m.B.copy_(B.to(DEVICE)) + m.W.copy_(base_keys(U, d // P).to(DEVICE)) + m.calibrate_power() + return m + + +def hadamard_book(vu=16, Lp=16): + B = torch.zeros(vu, Lp) + B[:, :vu] = torch.tensor(hadamard(vu).copy(), dtype=torch.float32) + return B + + +def random_ortho_book(vu=16, Lp=16, seed=7): + g = torch.Generator().manual_seed(seed) + A = torch.randn(Lp, Lp, generator=g) + Q, _ = torch.linalg.qr(A) + return Q[:vu].contiguous() + + +def report(name, m): + with torch.no_grad(): + Bn = m.unit_codebook() + G = Bn @ Bn.T + off = (G - torch.diag(torch.diag(G))).abs().max() + row = [name, "%.2e" % off] + for s in SNR: + row.append("%.4f" % ser(m, s, FRAMES)) + print("%-22s %-10s %-9s %-9s %-9s" % tuple(row)) + + +def main(): + print("%-22s %-10s %-9s %-9s %-9s" + % ("unit codebook", "max|off|", "0 dB", "10 dB", "20 dB")) + report("Walsh-Hadamard", fixed_model(hadamard_book())) + report("random orthogonal", fixed_model(random_ortho_book())) + print("%-22s %-10s %-9s %-9s %-9s" + % ("trained (paper)", "0.887", "0.8822", "0.2576", "0.0307")) + print("%-22s %-10s %-9s %-9s %-9s" + % ("OMA, resource matched", "-", + "%.4f" % L.oma_ser([0.0])[0], + "%.4f" % L.oma_ser([10.0])[0], + "%.4f" % L.oma_ser([20.0])[0])) + + + + + +def solo_check(): + """Splitting each candidate set from the superposition it must live + in. Orthogonal codewords are ideal for one user alone and are what + the single-user bound assumes, but the masked sets of different users + are then far from orthogonal to each other.""" + print() + print("%-22s %-12s %-12s" % ("unit codebook", "solo 10 dB", "4-user 10 dB")) + for name, B in (("random orthogonal", random_ortho_book()), + ("trained (retrain)", None)): + if B is None: + from exp_full import main_model + m = main_model() + else: + m = fixed_model(B) + print("%-22s %-12.4f %-12.4f" + % (name, ser(m, 10.0, FRAMES, solo=True), + ser(m, 10.0, FRAMES, solo=False))) + print("single-user ideal M-ary bound (separate MC): 0.1986") + + +if __name__ == "__main__": + main() + solo_check()