Revision history [back]
 | 1 | initial version |
Given the congruences:
(a*n-a+2) mod 83 = 0
(a*n-1884*a+119) mod 1683 = 0
(a*n-731*a+649) mod 5923 = 0
(a*n-6654*a+6572) mod 53307 = 0
(a*n-46652*a+46570) mod 373291 = 0
(a*n-46653*a+46571) mod 373299 = 0
and
|n|<|a|<27
ChatGPT wrote this for me: (Please excuse me for using AI)
# SageMath
# Reticolo bivariato LLL per le sei congruenze
# an + c_i*a + d_i == 0 (mod N_i)
#
# Bound: |n| < |a| < 27
#
# Nota: e' un metodo di tipo Coppersmith/Howgrave-Graham
# euristico per il sistema combinato, non una garanzia
# teorica di recupero delle radici.
from sage.all import *
# --------------------------------------------------
# 1. Dati
# --------------------------------------------------
Ns = [
ZZ(83),
ZZ(1683),
ZZ(5923),
ZZ(53307),
ZZ(373291),
ZZ(373299)
]
cs = [
ZZ(-1),
ZZ(-1884),
ZZ(-731),
ZZ(-6654),
ZZ(-46652),
ZZ(-46653)
]
ds = [
ZZ(2),
ZZ(119),
ZZ(649),
ZZ(6572),
ZZ(46570),
ZZ(46571)
]
X = ZZ(26)
Y = ZZ(25)
# Parametri del reticolo: aumentabili
t = 3
# --------------------------------------------------
# 2. Polinomi bivariati
# --------------------------------------------------
R.<a,n> = PolynomialRing(ZZ, order='lex')
f = [
a*n + cs[i]*a + ds[i]
for i in range(len(Ns))
]
print("Polinomi originali:")
for i in range(len(Ns)):
print("f_%d =" % i, f[i])
print("modulo", Ns[i])
# --------------------------------------------------
# 3. Combinazione modulare
# --------------------------------------------------
# P e' il prodotto dei moduli.
# Se tutte le congruenze sono soddisfatte, allora
# G(a,n) == 0 mod P.
#
# G = sum (P/N_i)*f_i.
#
# La condizione inversa non e' garantita:
# ogni candidato sara' verificato su tutte le f_i.
P = prod(Ns)
G = sum(
(P // Ns[i]) * f[i]
for i in range(len(Ns))
)
print("\nP =", P)
print("Grado totale di G =", G.total_degree())
# --------------------------------------------------
# 4. Generazione degli shift
# --------------------------------------------------
# Si usano shift monomiali di G e shift multipli di P.
# La base viene scalata con X^i * Y^j per incorporare
# i bound delle radici.
def exp_pairs(D):
return [
(i, j)
for total in range(D + 1)
for i in range(total + 1)
for j in [total - i]
]
shifts = []
for i, j in exp_pairs(t):
shifts.append(a**i * n**j * G)
for i, j in exp_pairs(t + 2):
shifts.append(P * a**i * n**j)
# Grado massimo dei polinomi di base
D = t + 4
mons_exp = exp_pairs(D)
# Vettori dei coefficienti scalati
rows = []
for h in shifts:
row = []
for i, j in mons_exp:
coeff = h.monomial_coefficient(a**i * n**j)
row.append(ZZ(coeff) * X**i * Y**j)
if any(v != 0 for v in row):
rows.append(row)
M = matrix(ZZ, rows)
print("\nDimensioni della matrice iniziale:", M.nrows(), "x", M.ncols())
# --------------------------------------------------
# 5. Riduzione LLL
# --------------------------------------------------
L = M.LLL()
print("Dimensioni della base ridotta:", L.nrows(), "x", L.ncols())
# Ricostruisce i polinomi dai vettori ridotti,
# rimuovendo la scalatura X^i*Y^j.
S.<u,v> = PolynomialRing(QQ, order='lex')
reduced_polys = []
for row in L.rows():
h = S.zero()
for k, (i, j) in enumerate(mons_exp):
scale = X**i * Y**j
if scale != 0 and row[k] != 0:
h += QQ(row[k]) / scale * u**i * v**j
if h != 0 and h not in reduced_polys:
reduced_polys.append(h)
print("\nPolinomi non nulli ottenuti:", len(reduced_polys))
for i, h in enumerate(reduced_polys[:8]):
print("\nh_%d =" % i, h)
# --------------------------------------------------
# 6. Risultanti per ricavare candidati per a
# --------------------------------------------------
T.<z> = PolynomialRing(QQ)
candidate_a = set()
# Cerca coppie di polinomi con risultante non nullo.
for i in range(len(reduced_polys)):
for j in range(i + 1, len(reduced_polys)):
h1 = reduced_polys[i]
h2 = reduced_polys[j]
try:
res = h1.resultant(h2, v)
except Exception:
continue
if res == 0:
continue
# Il risultante deve dipendere da u.
try:
rz = T(res(u=z))
except Exception:
continue
if rz == 0 or rz.degree() <= 0:
continue
# Le radici intere sono candidate per a.
try:
for root, mult in rz.roots(ZZ):
if abs(root) < 27 and root != 0:
candidate_a.add(ZZ(root))
except Exception:
pass
print("\nCandidati a dai risultanti:", sorted(candidate_a))
# --------------------------------------------------
# 7. Verifica rigorosa nelle congruenze originali
# --------------------------------------------------
def verifica(aa, nn):
if aa == 0:
return False
if not (abs(nn) < abs(aa) < 27):
return False
return all(
(aa*nn + cs[i]*aa + ds[i]) % Ns[i] == 0
for i in range(len(Ns))
)
soluzioni = set()
# Per ogni candidato a, n e' cercato solo nel bound
# ristretto. Questa fase filtra i candidati LLL.
for aa in candidate_a:
for nn in range(-25, 26):
if verifica(aa, ZZ(nn)):
soluzioni.add((aa, ZZ(nn)))
print("\nSoluzioni certificate trovate dal reticolo:")
for aa, nn in sorted(soluzioni):
print("a =", aa, ", n =", nn)
# --------------------------------------------------
# 8. Controllo indipendente entro i bound
# --------------------------------------------------
# Utile per distinguere un fallimento del reticolo
# da una vera assenza di soluzioni.
print("\nControllo esaustivo indipendente:")
soluzioni_esatte = []
for aa in range(-26, 27):
if aa == 0:
continue
for nn in range(-25, 26):
if abs(nn) < abs(aa) and verifica(ZZ(aa), ZZ(nn)):
soluzioni_esatte.append((ZZ(aa), ZZ(nn)))
print(soluzioni_esatte)
OUTPUT:
Polinomi originali:
f_0 = a*n - a + 2
modulo 83
f_1 = a*n - 1884*a + 119
modulo 1683
f_2 = a*n - 731*a + 649
modulo 5923
f_3 = a*n - 6654*a + 6572
modulo 53307
f_4 = a*n - 46652*a + 46570
modulo 373291
f_5 = a*n - 46653*a + 46571
modulo 373299
P = 6145999617209997401439561
Grado totale di G = 2
Dimensioni della matrice iniziale: 31 x 36
Dimensioni della base ridotta: 31 x 36
Polinomi non nulli ottenuti: 21
h_0 = -1279368*u*v - 11994075*u + 116902251
h_1 = 1279368*u*v^2 + 1759131*u*v - 95952600*u - 116902251*v + 935218008
h_2 = -1279368*u^2*v - 11994075*u^2 - 1279368*u*v + 104908176*u + 116902251
h_3 = 1279368*u*v^3 + 479763*u*v^2 - 14552811*u*v + 875567475*u - 116902251*v^2 + 1052120259*v - 8533864323
h_4 = 1279368*u^2*v^2 + 1759131*u^2*v - 95952600*u^2 - 28625859*u*v + 1762809183*u - 8066255319
h_5 = -1279368*u^3*v - 11994075*u^3 - 1279368*u^2*v + 104908176*u^2 - 1279368*u*v + 104908176*u + 116902251
h_6 = -1279368*u^2*v - 11994075*u^2 + 1279368*u*v^4 + 479763*u*v^3 - 4317867*u*v^2 + 132254667*u*v - 7751210949*u - 116902251*v^3 + 1052120259*v^2 - 9469082331*v + 76687876656
h_7 = -1279368*u^2*v^3 - 479763*u^2*v^2 + 18390915*u^2*v - 839585250*u^2 + 8155971*u*v^2 - 391806450*u*v + 23931378045*u + 9936691335*v - 153492655563
Candidati a dai risultanti: [-1, 17]
Soluzioni certificate trovate dal reticolo:
a = 17 , n = -4
Controllo esaustivo indipendente:
[(17, -4)]
This is exactly what I was looking for. Is this implementation correct? How could it be improved? What is the computational cost?
