def solve_quadratic_congruence(a, m):
"""
Находит все решения уравнения n^2 = a (mod m).
Args:
a (int): Остаток.
m (int): Модуль.
Returns:
list: Отсортированный список всех уникальных решений n в диапазоне [0, m-1].
"""
if m <= 1:
return []
a %= m
# Разложение m на простые множители вида p^k
prime_factorization = {}
d = 2
temp_m = m
while d * d <= temp_m:
while temp_m % d == 0:
prime_factorization[d] = prime_factorization.get(d, 0) + 1
temp_m //= d
d += 1
if temp_m > 1:
prime_factorization[temp_m] = prime_factorization.get(temp_m, 0) + 1
# Решения для каждого простого множителя p^k
solutions_by_modulus = {}
for p, k in prime_factorization.items():
pk = p ** k
solutions_by_modulus[pk] = _solve_for_prime_power(a % pk, p, k)
if not solutions_by_modulus[pk]:
# Если для одного из множителей нет решений, то и для m их нет
return []
# Комбинирование решений с помощью Китайской теоремы об остатках
# Инициализация с решениями для первого множителя
mods = list(solutions_by_modulus.keys())
if not mods:
return []
current_solutions = solutions_by_modulus[mods[0]]
current_mod = mods[0]
# Постепенное комбинирование решений
for i in range(1, len(mods)):
next_mod = mods[i]
next_sols = solutions_by_modulus[next_mod]
new_solutions = []
for sol1 in current_solutions:
for sol2 in next_sols:
# Решаем систему x = sol1 (mod current_mod), x = sol2 (mod next_mod)
res, mod_res = _crt_pair(sol1, current_mod, sol2, next_mod)
new_solutions.append(res)
current_solutions = new_solutions
current_mod *= next_mod
return sorted(list(set(current_solutions)))
def _solve_for_prime_power(a, p, k):
"""Решает n^2 = a (mod p^k)."""
if a == 0:
# Для n^2 = 0 (mod p^k), решения - это кратные p^ceil(k/2)
ceil_k_div_2 = (k + 1) // 2
root = pow(p, ceil_k_div_2)
return list(range(0, pow(p, k), root))
# Сначала решаем n^2 = a (mod p)
base_sols = _solve_for_prime(a % p, p)
if not base_sols:
return []
# Поднятие решений с mod p до mod p^k (Лемма Гензеля)
pk = pow(p, k)
all_sols = set()
for r in base_sols:
# Если p делит 2r, но не p^k, то нужен специальный подход (здесь для p=2)
if p == 2:
sols = _lift_hensel_p2(a, k, r)
all_sols.update(sols)
else: # Для нечетных простых
# f(x) = x^2 - a, f'(x) = 2x
# r_new = r - f(r) * (f'(r))^{-1} mod p^i
current_r = r
mod_p = p
for _ in range(1, k):
f_r = (current_r * current_r - a)
f_prime_r = 2 * current_r
# Находим (f'(r))^{-1} mod p
inv = _modInverse(f_prime_r, p)
# t = (f(r)/p^i) * (f'(r))^{-1} mod p
t = (f_r // mod_p) * inv
# r_{i+1} = r_i - t * p^i
current_r = (current_r - t * mod_p) % (mod_p * p)
mod_p *= p
all_sols.add(current_r)
return sorted(list(all_sols))
def _solve_for_prime(a, p):
"""Решает n^2 = a (mod p) с помощью алгоритма Тонелли-Шенкса."""
if p == 2:
return [a % 2]
# Проверка, является ли 'a' квадратичным вычетом по модулю p
if pow(a, (p - 1) // 2, p) != 1:
return [] if a % p != 0 else [0]
if a % p == 0:
return [0]
# Алгоритм Тонелли-Шенкса
if p % 4 == 3:
r = pow(a, (p + 1) // 4, p)
return [r, p - r]
# Разложение p-1 = Q * 2^S
Q = p - 1
S = 0
while Q % 2 == 0:
Q //= 2
S += 1
# Поиск неквадратичного вычета z
z = 2
while pow(z, (p - 1) // 2, p) == 1:
z += 1
M = S
c = pow(z, Q, p)
t = pow(a, Q, p)
R = pow(a, (Q + 1) // 2, p)
while t != 1:
if t == 0:
return [0]
# Поиск i, такого что t^(2^i) = 1 (mod p)
i = 0
temp_t = t
while temp_t != 1:
temp_t = (temp_t * temp_t) % p
i += 1
if i == M: # Не должно случиться, если 'a' - вычет
return []
b = pow(c, pow(2, M - i - 1, p - 1), p)
M = i
c = (b * b) % p
t = (t * c) % p
R = (R * b) % p
return [R, p - R]
def _lift_hensel_p2(a, k, r_base):
"""Поднятие решений для p=2."""
sols = {r_base}
mod = 2
for i in range(1, k):
new_sols = set()
mod_new = mod * 2
for r in sols:
if (r*r - a) % mod_new == 0:
new_sols.add(r)
if ((r+mod)*(r+mod) - a) % mod_new == 0:
new_sols.add(r + mod)
sols = new_sols
mod = mod_new
return list(sols)
def _extended_gcd(a, b):
"""Расширенный алгоритм Евклида."""
if a == 0:
return b, 0, 1
d, x1, y1 = _extended_gcd(b % a, a)
x = y1 - (b // a) * x1
y = x1
return d, x, y
def _modInverse(a, m):
"""Нахождение обратного по модулю."""
d, x, y = _extended_gcd(a, m)
if d != 1:
# Обратный элемент не существует
return None
return (x % m + m) % m
def _crt_pair(a1, m1, a2, m2):
"""Китайская теорема об остатках для двух уравнений."""
d, n1, n2 = _extended_gcd(m1, m2)
# x = a1 (mod m1), x = a2 (mod m2)
# Решение существует, если a1 = a2 (mod gcd(m1, m2))
if (a1 - a2) % d != 0:
return None, None # Нет решений
mod = m1 // d * m2
k = (a2 - a1) // d
x = a1 + k * n1 * m1
return x % mod, mod
# --- Примеры использования ---
if __name__ == '__main__':
# Пример 1: n^2 = 1 (mod 15)
# Решения: 1, 4, 11, 14
a1, m1 = 1, 15
solutions1 = solve_quadratic_congruence(a1, m1)
print(f"Решения для n^2 = {a1} (mod {m1}): {solutions1}")
# Пример 2: n^2 = 10 (mod 13)
# Решения: 6, 7
a2, m2 = 10, 13
solutions2 = solve_quadratic_congruence(a2, m2)
print(f"Решения для n^2 = {a2} (mod {m2}): {solutions2}")
# Пример 3: n^2 = 5 (mod 11)
# Решений нет
a3, m3 = 5, 11
solutions3 = solve_quadratic_congruence(a3, m3)
print(f"Решения для n^2 = {a3} (mod {m3}): {solutions3}")
# Пример 4: n^2 = 7 (mod 32)
# Решений нет, так как n^2 = 7 (mod 8) не имеет решений
a4, m4 = 7, 32
solutions4 = solve_quadratic_congruence(a4, m4)
print(f"Решения для n^2 = {a4} (mod {m4}): {solutions4}")
# Пример 5: n^2 = 1 (mod 8)
# Решения: 1, 3, 5, 7
a5, m5 = 1, 8
solutions5 = solve_quadratic_congruence(a5, m5)
print(f"Решения для n^2 = {a5} (mod {m5}): {solutions5}")
# Пример 6: n^2 = 71 (mod 105)
# Решения: 16, 19, 86, 89
a6, m6 = 71, 105
solutions6 = solve_quadratic_congruence(a6, m6)
print(f"Решения для n^2 = {a6} (mod {m6}): {solutions6}")
# Пример 7: n^2 = 0 (mod 100)
# Решения: 0, 10, 20, 30, 40, 50, 60, 70, 80, 90
a7, m7 = 0, 100
solutions7 = solve_quadratic_congruence(a7, m7)
print(f"Решения для n^2 = {a7} (mod {m7}): {solutions7}")