rorikstr/hyperion
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

sync

Browse source
This commit is contained in:
Rorik Star Platinum 2025-12-28 19:00:03 +03:00
parent 48f640b55d
commit 440a2fa01f
28 changed files with 1849 additions and 96 deletions
Download patch file Download diff file Expand all files Collapse all files

245
80-math/miller-rabin.md
View file

@ -2,30 +2,36 @@
пример кода, близкого к эталонному:
```rust
pub fn nth(n: u32) -> u32 {
let mut num = 1;
for _ in 0..=n {
loop {
num += 1;
if miller_rabin(num as u64) {
break;
}
}
}
num
fn main() {
let n: u64 = std::env::args().nth(1)
.and_then(|s| s.parse().ok())
.unwrap_or(35);
println!("{} -> {}", n, if is_prime(n) { "prime" } else { "composite" });
}
fn miller_rabin(n: u64) -> bool {
const HINT: &[u64] = &[2];
pub fn is_prime(n: u64) -> bool {
if n < 2 { return false; }
if n == 2 || n == 3 { return true; }
if n % 2 == 0 { return false; }
// we have a strict upper bound, so we can just use the witness
// table of Pomerance, Selfridge & Wagstaff and Jeaschke to be as
// efficient as possible, without having to fall back to
// randomness. Additional limits from Feitsma and Galway complete
// the entire range of `u64`. See also:
// https://en.wikipedia.org/wiki/Miller%E2%80%93Rabin_primality_test#Testing_against_small_sets_of_bases
let s = (n - 1).trailing_zeros();
let d = (n - 1) >> s;
let witnesses = get_witnesses(n);
for &a in witnesses {
if !miller_rabin_test(a, d, n, s) {
return false;
}
}
true
}
// O(1) lookup based on input size
fn get_witnesses(n: u64) -> &'static [u64] {
const WITNESSES: &[(u64, &[u64])] = &[
(2_046, HINT),
(2_046, &[2]),
(1_373_652, &[2, 3]),
(9_080_190, &[31, 73]),
(25_326_000, &[2, 3, 5]),
@ -35,88 +41,50 @@ fn miller_rabin(n: u64) -> bool {
(3_474_749_660_382, &[2, 3, 5, 7, 11, 13]),
(341_550_071_728_320, &[2, 3, 5, 7, 11, 13, 17]),
(3_825_123_056_546_413_050, &[2, 3, 5, 7, 11, 13, 17, 19, 23]),
(std::u64::MAX, &[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37]),
(u64::MAX, &[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37]),
];
if n % 2 == 0 {
return n == 2;
}
if n == 1 {
return false;
}
let mut d = n - 1;
let mut s = 0;
while d % 2 == 0 {
d /= 2;
s += 1
}
let witnesses = WITNESSES
.iter()
.find(|&&(hi, _)| hi >= n)
.map(|&(_, wtnss)| wtnss)
.unwrap();
'next_witness: for &a in witnesses.iter() {
let mut power = mod_exp(a, d, n);
assert!(power < n);
if power == 1 || power == n - 1 {
continue 'next_witness;
}
for _r in 0..s {
power = mod_sqr(power, n);
assert!(power < n);
if power == 1 {
return false;
}
if power == n - 1 {
continue 'next_witness;
}
}
return false;
}
true
WITNESSES.iter()
.find(|(limit, _)| *limit >= n)
.map(|(_, w)| *w)
.unwrap()
}
fn mod_mul_(a: u64, b: u64, m: u64) -> u64 {
(u128::from(a) * u128::from(b) % u128::from(m)) as u64
fn miller_rabin_test(a: u64, d: u64, n: u64, s: u32) -> bool {
let mut x = mod_pow(a, d, n);
if x == 1 || x == n - 1 {
return true;
}
for _ in 1..s {
x = mod_mul(x, x, n);
if x == n - 1 {
return true;
}
}
false
}
// === Math Core ===
#[inline]
fn mod_mul(a: u64, b: u64, m: u64) -> u64 {
match a.checked_mul(b) {
Some(r) => {
if r >= m {
r % m
} else {
r
}
}
None => mod_mul_(a, b, m),
}
((a as u128 * b as u128) % m as u128) as u64
}
fn mod_sqr(a: u64, m: u64) -> u64 {
if a < (1 << 32) {
let r = a * a;
if r >= m { r % m } else { r }
} else {
mod_mul_(a, a, m)
fn mod_pow(mut base: u64, mut exp: u64, m: u64) -> u64 {
let mut res = 1;
base %= m;
while exp > 0 {
if exp % 2 == 1 {
res = mod_mul(res, base, m);
}
base = mod_mul(base, base, m);
exp /= 2;
}
res
}
fn mod_exp(mut x: u64, mut d: u64, n: u64) -> u64 {
let mut ret: u64 = 1;
while d != 0 {
if d % 2 == 1 {
ret = mod_mul(ret, x, n)
}
d /= 2;
x = mod_sqr(x, n);
}
ret
}
```
@ -153,7 +121,9 @@ let mut s = 0;
while d % 2 == 0 { d /= 2; s += 1 }
```
**Математика:** Каждое простое удовлетворяет **малой теореме Ферма**: \(a^{n-1} \equiv 1 \pmod n\).
**Математика:** Каждое простое удовлетворяет **малой теореме Ферма**:
\(a^{n-1} \equiv 1 \pmod n\).
Но вместо прямого вычисления \(a^{n-1}\) мы раскладываем показатель:
\[ n-1 = d \cdot 2^s \]
@ -360,3 +330,96 @@ fn fermat_test(n: u64, a: u64) -> bool {
[8](https://journalpro.ru/articles/malaya-teorema-ferma-i-ee-primenenie-v-kriptosistemakh/)
[9](https://www.reddit.com/r/math/comments/2v6db7/intuition_for_fermats_little_theorem/)
[10](https://math.fandom.com/ru/wiki/%D0%9C%D0%B0%D0%BB%D0%B0%D1%8F_%D1%82%D0%B5%D0%BE%D1%80%D0%B5%D0%BC%D0%B0_%D0%A4%D0%B5%D1%80%D0%BC%D0%B0)
Ты нащупал самую суть! Это называется **"Binary Exponentiation" (Бинарное возведение в степень)**.
Твоя догадка: _"мы ищем минимальную степень... домножаем... переиспользуем"_ — абсолютно верна.
Давай разложим это на механику "Сборщика" и "Генератора".
## 1. Генератор: Как мы переиспользуем `base`?
Посмотри на строку:
rust
`base = mod_mul(base, base, m);`
В каждой итерации цикла эта строка делает одно и то же: **возводит текущую базу в квадрат**.
Представь, что на входе `base` — это число $X$.
- **Итерация 1:** `base` становится $X^2$.
- **Итерация 2:** Мы берем этот $X^2$ и умножаем сам на себя: $(X^2) \cdot (X^2) = X^4$.
- **Итерация 3:** Берем $X^4$ и умножаем сам на себя: $(X^4) \cdot (X^4) = X^8$.
- **Итерация 4:** ... $X^{16}$.
**Вот оно, переиспользование!** Вместо того чтобы умножать $X \cdot X \cdot X \dots$ (линейно), мы прыгаем по степеням двойки: 1, 2, 4, 8, 16, 32...
Мы **генерируем** эти "кирпичики" ($X^1, X^2, X^4, X^8$) на лету.
## 2. Сборщик: Как мы собираем результат `res`?
Теперь посмотри на условие:
rust
`if exp % 2 == 1 { res = mod_mul(res, base, m); }`
Здесь мы решаем: **нужен нам этот кирпичик или нет?**
Допустим, нам нужно вычислить $3^{13}$.
Двоичное представление 13 — это `1101`.
Это значит:
13=1⋅8+1⋅4+0⋅2+1⋅113 = 1 \cdot 8 + 1 \cdot 4 + 0 \cdot 2 + 1 \cdot 113=1⋅8+1⋅4+0⋅2+1⋅1
То есть:
313=38⋅34⋅(пропускаем 32)⋅313^{13} = 3^8 \cdot 3^4 \cdot (\text{пропускаем } 3^2) \cdot 3^1313=38⋅34⋅(пропускаем 32)⋅31
Мы просто берем только те степени, где в двоичном коде стоит **1**.
## Визуализация процесса (Трассировка)
Давай пройдем этот цикл для $3^{13} \pmod{100}$.
Изначально: `res = 1`, `base = 3`. `exp = 13` (`1101` bin).
| Итерация | exp (dec) | exp (bin) | Бит (четность) | Действие с `res` (Сборщик) | Состояние `res` | Действие с `base` (Генератор) | Новое `base` |
| -------- | --------- | --------- | -------------- | -------------------------- | ---------------------------------- | ----------------------------- | ------------------- |
| **1** | 13 | `...1101` | **1 (нечет)** | Нужен! Умножаем на $3^1$ | $1 \cdot 3 = \mathbf{3}$ | $3 \cdot 3$ (квадратируем) | $3^2 = 9$ |
| **2** | 6 | `...0110` | **0 (чет)** | Не нужен. Пропускаем. | $\mathbf{3}$ | $9 \cdot 9$ (квадратируем) | $3^4 = 81$ |
| **3** | 3 | `...0011` | **1 (нечет)** | Нужен! Умножаем на $3^4$ | $3 \cdot 81 = 243 \to \mathbf{43}$ | $81 \cdot 81 \to 61$ | $3^8 = 61$ |
| **4** | 1 | `...0001` | **1 (нечет)** | Нужен! Умножаем на $3^8$ | $43 \cdot 61 \to \mathbf{23}$ | $61 \cdot 61 \dots$ | $3^{16}$ (не важно) |
**Итог:** Мы получили 23.
Проверка: $3^{13} = 1,594,323$.
$1,594,323 \pmod{100} = 23$. **Сошлось!**
## Доказательство формулы остатков
Ты спрашиваешь: _"Как мы доказываем, что взятие остатка на каждом шаге не ломает математику?"_
Это фундаментальное свойство модульной арифметики:
**Остаток произведения равен произведению остатков.**
(A⋅B)(modM)=((A(modM))⋅(B(modM)))(modM)(A \cdot B) \pmod M = ((A \pmod M) \cdot (B \pmod M)) \pmod M(A⋅B)(modM)=((A(modM))⋅(B(modM)))(modM)
**Пример на пальцах:**
Пусть $M=10$ (мы смотрим только на последнюю цифру).
Мы хотим узнать последнюю цифру $12 \cdot 13$.
1. **По-честному:** $12 \cdot 13 = 156$. Последняя цифра **6**.
2. **По-модульному:**
- От 12 берем остаток $\to 2$.
- От 13 берем остаток $\to 3$.
- Умножаем остатки: $2 \cdot 3 = 6$.
- Результат тот же: **6**.
Именно это правило разрешает нам делать `base %= m` и `mod_mul` на каждом шаге. Мы просто отбрасываем "лишние круги", которые никак не влияют на конечный остаток. Без этого правила вся современная криптография (RSA, Diffie-Hellman) просто не работала бы на компьютерах.