def prime_factors(n, return_e: bool = True): """ nの素因数分解を指数付きで返す 計算量O(√N) O(n^(1/4))もある https://qiita.com/drken/items/a14e9af0ca2d857dad23#4-%E7%B4%A0%E5%9B%A0%E6%95%B0%E5%88%86%E8%A7%A3 """ assert 1 <= n pfs = [] a = 2 while a**2 <= n: if n % a == 0: ex = 0 while n % a == 0: ex += 1 n //= a pfs.append((a, ex) if return_e else a) a += 1 if n != 1: pfs.append((n, 1) if return_e else n) return pfs def primitive_root(p: int) -> int: """素数pに対する何らかの原始根rを一つ求める (rの位数はp-1) https://manabitimes.jp/math/842 https://37zigen.com/primitive-root/ 最速ではない """ if p == 2: return 1 from random import shuffle p_1_pfs = prime_factors(p-1, return_e=False) candidates_for_r = list(range(2, p)) shuffle(candidates_for_r) # 乱択で原始根を探す for cfr in candidates_for_r: for pf in p_1_pfs: if pow(cfr, (p-1)//pf, p) == 1: # p-1//pfが位数候補 break else: # 全ての位数候補に対して、その否定がなされた。 return cfr else: raise AssertionError def discrete_logarithm(x, y, mod): """離散対数問題 x^k≡yとなるようなkを求める O(√mod)""" babys = {} sqrt = int(mod**0.5)+1 # baby-step x_baby = 1 for i in range(sqrt): babys[x_baby] = i # x^(sqrt-1)まで求めておく if x_baby == y: return i x_baby = x_baby * x % mod # giant-step x_giant = pow(x_baby, mod-2, mod) # x^(-sqrt) for i in range(1, sqrt+1): y = y * x_giant % mod if y in babys: # ここで√modずつ見ていく為に速い return babys[y] + i*sqrt return -1 # 原始根でない限り、必ず存在するとは言えない """ 関連する話題 これはメモ20210505 フェルマーの小定理 pが素数 aがpの倍数でない正の整数の時 a^(p-1)≡1 (mod p) 証明 a^p=aを示す 帰納法にすると(1+a)^p≡1+a^p≡1+a (∵二項係数) よって示された 別解 ia(i={1,2,…,p-1})にはmod pで1,2,…,p-1が出る(∵背理法) よってa^(p-1)*(p-1)!≡(p-1)!⇔a^(p-1)≡1(∵(p-1)!とpは互いに素) https://manabitimes.jp/math/680 https://qiita.com/drken/items/6b4031ccbb2cab7436f3 レプユニット数(11…11) これは完全に数学 https://manabitimes.jp/math/1374 NTT(整数環FFT)(Number Theoretic Transform) FFTで出てくる丸め誤差の問題を解決する畳み込みの手法 https://math314.hateblo.jp/entry/2015/05/07/014908 garnerのアルゴリズム 本質的にはあまり変わらない https://www.csee.umbc.edu/~lomonaco/s08/441/handouts/Garner-Alg-Example.pdf """ def _inv_gcd(a, b) -> tuple: """ 返り値は(gcd(a,b),x) (ただしxa≡g (mod b) 0<=x=0) m0 += b // s return (s, m0) def eea(a: int, b: int) -> tuple: """ ax+by=gcd(a,b)なる(x,y)を求める(拡張されたユークリッドの互除法) (extended_euclidean_algorithm)(元のユークリッドの互除法はgcdを求める手法を指す) https://qiita.com/drken/items/b97ff231e43bce50199a https://ja.wikipedia.org/wiki/ユークリッドの互除法 <-英語版が優秀すぎ https://en.wikipedia.org/wiki/Extended_Euclidean_algorithm また、以上の話題は連分数展開とも関連があるらしい """ g, x = _inv_gcd(a, abs(b)) y = (b//abs(b))*(g-a*x)//abs(b) assert (g-a*x) % abs(b) == 0 return (x, y) def crt(rs: list, mods: list) -> tuple: """ 同じ長さのリストr,mを引数に取り、このリストの長さをnとした時 x≡rs[i] (mod ms[i]) ∀i∊{0,1,…,n-1}を解く 答えが存在するならばx≡y(mod z) (0<=y