# 从辗转相除法到求逆元，数论算法初体验

## 暴力解法

def gcd(a, b):
ret = 0
for i in range(min(a, b)):
if a % i == 0 and b % i == 0:
ret = i
return ret

## 辗转相除法

$gcd(a, b) = gcd(b, r), a = bq + r$

\begin{aligned} a = su, b = tu \\ r = a – bq = su – tuq = (s – tq) u\\ \end{aligned}

\begin{aligned} b = sv, r = tv\\ a = bq + r = svq + tv = v(sq + t) \end{aligned}

def gcd(a, b):
if a < b:
a, b = b, a

if a % b == 0:
return b
return gcd(b, a % b)

def gcd(a, b):
return a if b == 0 else gcd(b, a % b)

## 拓展欧几里得

\begin{aligned} a*(x_0 + (b / gcd) * t) + b*(yo-(a/gcd)*t) \\ a*x_0+ b*y_0 + abt / gcd – abt/gcd = gcd \end{aligned}

\begin{aligned} \end{aligned}

\begin{aligned} gcd &= b*x_0 + (a\%b)*y_0 \\ &= b*x_0 + (a – (a/b)*b)*y_0 \\ &= b*x_0 + a*y_0 – (a/b)*b*y_0 \\ &= a*y_0 + b*(x_0 – (a/b)*b*y_0) \end{aligned}

def exgcd(a, b, x=1, y=0):
# 当b=0的时候return
if b == 0:
return a, x, y
# 递归调用，获取b, a%b时的gcd与通项解
gcd, x, y = exgcd(b, a%b, x, y)
# 代入，得到新的通项解
x, y = y, x - a//b*y
return gcd, x, y

## 逆元与解逆元

\begin{aligned} ax + by &= 1\\ ax \% b + by \% b &= 1 \% b\\ ax\%b &= 1\%b\\ ax &= 1 \pmod b \end{aligned}

def cal_inv(a, m):
gcd, x, y = exgcd(a, m)
# 如果gcd不为1，那么说明没有逆元，返回-1
return (x % m + m) % m if gcd == 1 else -1

$a^{m-1}\equiv 1 \pmod m$

$a^{m-2} \equiv inv(a) \pmod m$