oil_spill/recover_o.sage.py
2023-07-02 03:26:39 -04:00

165 lines
4.4 KiB
Python

# This file was *autogenerated* from the file recover_o.sage
from sage.all_cmdline import * # import sage library
_sage_const_3301 = Integer(3301); _sage_const_44 = Integer(44); _sage_const_112 = Integer(112); _sage_const_2 = Integer(2); _sage_const_1 = Integer(1); _sage_const_4 = Integer(4); _sage_const_0 = Integer(0); _sage_const_10 = Integer(10); _sage_const_200 = Integer(200)
from random import randint
from tqdm import tqdm
F = GF(_sage_const_3301 , names=("z", ))
f_a = list(F)
# m = 32
# n = 80
m = _sage_const_44
n = _sage_const_112
# m = 72
# n = 184
# m = 8
# n = 20
# m = 2
# n = 5
O = random_matrix(F, (n - m), m)
output = ""
poly_m = []
d = n / (m * log(_sage_const_3301 , _sage_const_2 ).n())
N = ceil(sqrt((n + _sage_const_1 )/_sage_const_4 )) + _sage_const_1
print("d: ", d)
print("N: ", N)
z = zero_matrix(F, m, (n - m))
for i in tqdm(range(m)):
P1 = random_matrix(F, (n - m), (n - m))
for j in range(_sage_const_0 , len(P1.rows())):
for k in range(_sage_const_0 , j):
P1[j, k] = _sage_const_0
P2 = random_matrix(F, (n - m), m)
P3 = -O.T * P1 * O - O.T * P2
for j in range(_sage_const_0 , len(P3.rows())):
for k in range(j+_sage_const_1 , len(P3.rows())):
P3[j, k] += P3[k, j]
P3[k, j] = _sage_const_0
for i in P1:
for j in i:
output += hex(f_a.index(j))[_sage_const_2 :]
for i in P2:
for j in i:
output += hex(f_a.index(j))[_sage_const_2 :]
for i in P3:
for j in i:
output += hex(f_a.index(j))[_sage_const_2 :]
P = block_matrix([ [P1, P2], [z, P3]])
poly_m.append(P)
v = matrix(F, n, _sage_const_1 , [randint(_sage_const_0 , _sage_const_1 ) for i in range(n)])
oil_basis = block_matrix(F, _sage_const_2 , _sage_const_1 , [O, identity_matrix(F, m)])
o1 = zero_matrix(F, n, _sage_const_1 )
o2 = zero_matrix(F, n, _sage_const_1 )
for i in oil_basis.columns():
o1 += F.random_element() * matrix(F, n, _sage_const_1 , list(i))
o2 += F.random_element() * matrix(F, n, _sage_const_1 , list(i))
o_hat1 = v + o1
# o_hat2 = v + o2
print("o1 :", [i[_sage_const_0 ] for i in o1])
# print("o2 :", [i[0] for i in o2])
print("v :", [i[_sage_const_0 ] for i in v])
print("o1+v :", [i[_sage_const_0 ] for i in o_hat1])
# print("o2+v :", [i[0] for i in o_hat2])
print()
# for P in poly_m:
# print("oPo:", o1.T * P * o1)
# print("vPv:", v.T * P * v)
# print("qPq:", o_hat.T * P * o_hat)
# correction = o_hat.T * (P.T + P) * v
# print("oPx:", correction)
# print("tot:", v.T * P * v + o_hat.T * P * o_hat - correction)
# print()
# print()
# P = poly_m[0]
vec = [[(o_hat1.T * (P.T + P))[_sage_const_0 ][i] for i in range(n)] for P in poly_m]
s = [(v.T * P * v + o_hat1.T * P * o_hat1)[_sage_const_0 ][_sage_const_0 ] for P in poly_m]
# print("s: ", (v.T * P * v + o_hat.T * P * o_hat)[0][0])
# print("k: ", vec)
# out = 0
# for i in range(len(vec)):
# out += int(vec[i]) * int(v[i][0])
#
# t = matrix(QQ, 1, len(vec) + 2, [int(i[0]) for i in v] + [-(out // 3301), -1])
# print("t: ", t)
d = n / (m * log(_sage_const_3301 , _sage_const_2 ).n())
N = ceil(sqrt((n + _sage_const_1 )/_sage_const_4 )) + _sage_const_1
print("d: ", d)
print("N: ", N)
b1 = identity_matrix(QQ, n)
b2 = zero_matrix(QQ, n, _sage_const_1 )
b3 = matrix(QQ, m, n, [[N * int(i) for i in j] for j in vec]).T
b4 = zero_matrix(QQ, m, n)
b5 = zero_matrix(QQ, m, _sage_const_1 )
b6 = N * _sage_const_3301 * identity_matrix(QQ, m)
b7 = matrix(QQ, _sage_const_1 , n, [_sage_const_1 /_sage_const_2 ]*n)
b8 = matrix(QQ, _sage_const_1 , _sage_const_1 , [_sage_const_1 /_sage_const_2 ])
b9 = matrix(QQ, _sage_const_1 , m, [N * int(i) for i in s])
B = block_matrix([[b1, b2, b3], [b4, b5, b6], [b7, b8, b9]])
print(B)
print()
BH = B.LLL()
BH = matrix(ZZ, _sage_const_2 * BH)
print(BH)
for beta in tqdm(range(_sage_const_0 , _sage_const_10 )):
BH = BH.BKZ(block_size=beta + _sage_const_1 , fp='rr', precision=_sage_const_200 )
BH = _sage_const_1 /_sage_const_2 * matrix(QQ, BH)
print(BH)
print()
for out_vec in BH:
if abs(out_vec[-(m + _sage_const_1 )]) == _sage_const_1 /_sage_const_2 :
if out_vec[-(m + _sage_const_1 )] == -_sage_const_1 /_sage_const_2 :
testvec = [i + _sage_const_1 /_sage_const_2 for i in out_vec]
else:
testvec = [i + _sage_const_1 /_sage_const_2 for i in -out_vec]
print(testvec)
print()
print([i[_sage_const_0 ] for i in v])