fn solve_linear_congruence(a: u64, b: u64, m: u64) -> Option<(u64, u64)>Examples found in repository?
crates/competitive/src/math/discrete_logarithm.rs (line 333)
313 fn discrete_logarithm(
314 &self,
315 a: u64,
316 b: u64,
317 br_primes: &[BarrettReduction<u64>],
318 ) -> Option<(u64, u64)> {
319 let p = self.p;
320 let ord = self.ord;
321 let br = BarrettReduction::<u128>::new(p as u128);
322 let a = br.rem(a as _) as u64;
323 let b = br.rem(b as _) as u64;
324 if a == 0 {
325 return if b == 0 { Some((1, 1)) } else { None };
326 }
327 if b == 0 {
328 return None;
329 }
330
331 let x = self.index_calculus(a, br_primes)?;
332 let y = self.index_calculus(b, br_primes)?;
333 solve_linear_congruence(x, y, ord)
334 }
335}
336
337thread_local!(
338 static IC: UnsafeCell<IndexCalculus> = UnsafeCell::new(IndexCalculus::new());
339);
340
341pub fn discrete_logarithm_prime_mod(a: u64, b: u64, p: u64) -> Option<u64> {
342 IC.with(|ic| unsafe { &mut *ic.get() }.discrete_logarithm(a, b, p))
343 .map(|t| t.0)
344}
345
346/// a^x ≡ b (mod n), a has order p^e
347fn pohlig_hellman_prime_power_order(a: u64, b: u64, n: u64, p: u64, e: u32) -> Option<u64> {
348 let br = BarrettReduction::<u128>::new(n as u128);
349 let mul = |x: u64, y: u64| br.rem(x as u128 * y as u128) as u64;
350 let block_size = (p as f64).sqrt().ceil() as u64;
351 let mut baby = HashMap::<u64, u64>::new();
352 let g = pow(a, p.pow(e - 1), &br);
353 let mut xj = 1;
354 for j in 0..block_size {
355 baby.entry(xj).or_insert(j);
356 xj = mul(xj, g);
357 }
358 let xi = modinv(xj, n);
359 let mut t = 0u64;
360 for k in 0..e {
361 let mut h = pow(mul(modinv(pow(a, t, &br), n), b), p.pow(e - 1 - k), &br);
362 let mut ok = false;
363 for i in (0..block_size * block_size).step_by(block_size as usize) {
364 if let Some(j) = baby.get(&h) {
365 t += (i + j) * p.pow(k);
366 ok = true;
367 break;
368 }
369 h = mul(h, xi);
370 }
371 if !ok {
372 return None;
373 }
374 }
375 Some(t)
376}
377
378/// a^x ≡ b (mod p^e)
379fn discrete_logarithm_prime_power(a: u64, b: u64, p: u64, e: u32) -> Option<(u64, u64)> {
380 assert_ne!(p, 0);
381 assert_ne!(e, 0);
382 let n = p.pow(e);
383 assert!(a < n);
384 assert!(b < n);
385 assert_eq!(gcd(a, p), 1);
386 if p == 1 {
387 return Some((0, 1));
388 }
389 if a == 0 {
390 return if b == 0 { Some((1, 1)) } else { None };
391 }
392 if b == 0 {
393 return None;
394 }
395 if e == 1 {
396 return IC.with(|ic| unsafe { &mut *ic.get() }.discrete_logarithm(a, b, p));
397 }
398 let br = BarrettReduction::<u128>::new(n as _);
399 if p == 2 {
400 if e >= 3 {
401 if a % 4 == 1 && b % 4 != 1 {
402 return None;
403 }
404 let aa = if a % 4 == 1 { a } else { n - a };
405 let bb = if b % 4 == 1 { b } else { n - b };
406 let g = 5;
407 let ord = n / 4;
408 let x = pohlig_hellman_prime_power_order(g, aa, n, p, e - 2)?;
409 let y = pohlig_hellman_prime_power_order(g, bb, n, p, e - 2)?;
410 let t = solve_linear_congruence(x, y, ord)?;
411 match (a % 4 == 1, b % 4 == 1) {
412 (true, true) => Some(t),
413 (false, true) if t.0 % 2 == 0 => Some((t.0, lcm(t.1, 2))),
414 (false, false) if t.0 % 2 == 1 => Some((t.0, lcm(t.1, 2))),
415 (false, false) if a == b => Some((1, lcm(t.1, 2))),
416 _ => None,
417 }
418 } else if a == 1 {
419 if b == 1 { Some((0, 1)) } else { None }
420 } else {
421 assert_eq!(a, 3);
422 if b == 1 {
423 Some((0, 2))
424 } else if b == 3 {
425 Some((1, 2))
426 } else {
427 None
428 }
429 }
430 } else {
431 let ord = n - n / p;
432 let pf_ord = prime_factors(ord);
433 let g = (2..)
434 .find(|&g| check_primitive_root(g, ord, &br, &pf_ord))
435 .unwrap();
436 let mut pf_p = prime_factors(p - 1);
437 pf_p.push((p, e - 1));
438 let mut abm = vec![];
439 for (q, c) in pf_p {
440 let m = q.pow(c);
441 let d = ord / m;
442 let gg = pow(g, d, &br);
443 let aa = pow(a, d, &br);
444 let bb = pow(b, d, &br);
445 let x = pohlig_hellman_prime_power_order(gg, aa, n, q, c)?;
446 let y = pohlig_hellman_prime_power_order(gg, bb, n, q, c)?;
447 abm.push((x, y, m));
448 }
449 solve_linear_congruences(abm)
450 }
451}