Streaming Polyphase Factor Decomposition / SKAT

Download original file

This is the original source text. Reading it does not run or verify the program.

Read the typeset paper

1 import unittest
2 import numpy as np
3 
4 from approximate_lightcone import (
5     apply_polynomial,
6     block_transfer_matrix,
7     phase_norm_vector,
8     polynomial_difference_envelope,
9     power_difference_envelope,
10     power_phase_envelope,
11     threshold_operator,
12 )
13 
14 
15 class ApproximateLightConeTests(unittest.TestCase):
16     @classmethod
17     def setUpClass(cls):
18         rng = np.random.default_rng(20260923)
19         cls.q = 9
20         cls.r = 4
21         cls.widths = [cls.r] * cls.q
22         n = cls.q * cls.r
23 
24         L = np.zeros((n, n))
25         for h in range(cls.q):
26             sh = slice(h * cls.r, (h + 1) * cls.r)
27             L[sh, sh] = 2.1 * np.eye(cls.r) + 0.04 * rng.normal(size=(cls.r, cls.r))
28             if h + 1 < cls.q:
29                 sg = slice((h + 1) * cls.r, (h + 2) * cls.r)
30                 L[sh, sg] = 0.10 * rng.normal(size=(cls.r, cls.r))
31         H0 = L.T @ L
32 
33         E = np.zeros_like(H0)
34         for h in range(cls.q):
35             for g in range(h + 2, cls.q):
36                 sh = slice(h * cls.r, (h + 1) * cls.r)
37                 sg = slice(g * cls.r, (g + 1) * cls.r)
38                 C = rng.normal(size=(cls.r, cls.r))
39                 E[sh, sg] = C
40                 E[sg, sh] = C.T
41         E /= np.linalg.norm(E, ord=2)
42         H = H0 + 0.04 * E
43         scale = np.linalg.norm(H, ord=2)
44         cls.H = H / scale
45         cls.H0 = H0 / scale
46         cls.B = block_transfer_matrix(cls.H, cls.widths)
47         cls.B0 = block_transfer_matrix(cls.H0, cls.widths)
48         cls.X = np.zeros((n, 3))
49         cls.X[:cls.r] = rng.normal(size=(cls.r, 3))
50         cls.rho = phase_norm_vector(cls.X, cls.widths)
51 
52     def test_phase_envelope(self):
53         Y = self.X.copy()
54         for k in range(1, 8):
55             Y = self.H @ Y
56             actual = phase_norm_vector(Y, self.widths)
57             bound = power_phase_envelope(self.B, self.rho, k)
58             self.assertTrue(np.all(actual <= bound + 1e-12))
59 
60     def test_thresholded_power_difference(self):
61         Y = self.X.copy()
62         Y0 = self.X.copy()
63         for k in range(1, 8):
64             Y = self.H @ Y
65             Y0 = self.H0 @ Y0
66             actual = phase_norm_vector(Y - Y0, self.widths)
67             bound = power_difference_envelope(self.B, self.B0, self.rho, k)
68             self.assertTrue(np.all(actual <= bound + 1e-12))
69 
70     def test_polynomial_difference(self):
71         coeff = [0.2, -0.8, 1.3, -0.7, 0.25]
72         actual = phase_norm_vector(
73             apply_polynomial(self.H, self.X, coeff)
74             - apply_polynomial(self.H0, self.X, coeff),
75             self.widths,
76         )
77         bound = polynomial_difference_envelope(
78             self.B, self.B0, self.rho, coeff
79         )
80         self.assertTrue(np.all(actual <= bound + 1e-12))
81 
82     def test_hard_threshold_bound(self):
83         Ht, Bt, _ = threshold_operator(self.H, self.widths, 0.003)
84         k = 5
85         actual = phase_norm_vector(
86             np.linalg.matrix_power(self.H, k) @ self.X
87             - np.linalg.matrix_power(Ht, k) @ self.X,
88             self.widths,
89         )
90         bound = power_difference_envelope(self.B, Bt, self.rho, k)
91         self.assertTrue(np.all(actual <= bound + 1e-12))
92 
93     def test_exact_case_recovers_zero_leakage(self):
94         # H0 is block tridiagonal and X starts in phase zero.
95         Y = self.X.copy()
96         for k in range(1, 6):
97             Y = self.H0 @ Y
98             norms = phase_norm_vector(Y, self.widths)
99             self.assertLessEqual(np.max(norms[k + 1:]), 1e-13)
100 
101 
102 if __name__ == "__main__":
103     unittest.main(verbosity=2)

Back to source controls

Back to Streaming Polyphase Factor Decomposition / SKAT