拓展欧几里德算法

拓展欧几里德算法

定义

已知整数a，b，必存在整数x，y，使得
ax + by = gcd(a,b).

例子

用类似辗转相除法，求二元一次不定方程 47x+30y=1的整数解。

• 47 = 30 $\times$ 1 + 17
• 30 = 17 $\times$ 1 + 13
• 17 = 13 $\times$ 1 + 4
• 13 = 4 $\times$ 3 + 1

然后把它们改写成“余数等于”的形式

• 17 = 47 $\times$ 1 + 30 $\times$ (-1) //式1
• 13 = 30 $\times$ 1 + 17 $\times$ (-1) //式2
• 4 = 17 $\times$ 1 + 13 $\times$ (-1) //式3
• 1 = 13 $\times$ 1 + 4 $\times$ (-3)

然后把它们“倒回去”

• 1 = 13 $\times$ 1 + 4 $\times$ (-3)
• 1 = 13 $\times 1$ + [17 $\times$ 1 + 13 $\times$ (-1)] * (-3) //应用式3
• 1 = 17 $\times$ (-3) + 13 $\times$ 4
• 1 = 17 $\times$ (-3) + [30 $\times$ 1 + 17 $\times$ (-1)] $\times$ 4 //应用式2
• 1 = 30 $\times$ 4 + 17 $\times$ (-7)
• 1 = 30 $\times$ 4 + [47 $\times$ 1 + 30 $\times$ (-1)] $\times$ (-7) //应用式1
• 1 = 47 $\times$ (-7) + 30 $\times$ 11
  得解 x=-7, y=11 。

这个过程可以用矩阵表示（其中q表示商，r表示余数）

$\begin{pmatrix} a \ b \end{pmatrix} = \prod_{i=0}^{N} \begin{pmatrix} q_{i} & 1 \\ 1 & 0 \end{pmatrix} \begin{pmatrix} r_{N-1} \\ 0 \end{pmatrix}$

$\begin{pmatrix} 47 \ 30 \end{pmatrix} = \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}\begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix} \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}\begin{pmatrix} 3 & 1 \\ 1 & 0 \end{pmatrix} \begin{pmatrix} 4 & 1 \\ 1 & 0 \end{pmatrix} \begin{pmatrix} 1 \\ 0 \end{pmatrix} =$$\begin{pmatrix} 47 & 11 \\ 30 & 7 \end{pmatrix} \begin{pmatrix} 1 \\ 0 \end{pmatrix}\begin{pmatrix}1 \\ 0 \end{pmatrix}= \begin{pmatrix} -7 & 11 \\ 30 & -47 \end{pmatrix} \begin{pmatrix} 47 \\ 30 \end{pmatrix}$

或者用[[初等变换]]
$\begin{pmatrix}47 & 30\\1 & 0\\0 & 1\end{pmatrix} \rightarrow \begin{pmatrix}17 & 30\\1 & 0\\-1 & 1\end{pmatrix} \rightarrow \begin{pmatrix}17 & 13\\1 & -1\\-1 & 2\end{pmatrix}$
$\rightarrow \begin{pmatrix}4 & 13\\2 & -1\\-3 & 2\end{pmatrix} \rightarrow \begin{pmatrix}4 & 1\\2 & -7\\-3 & 11\end{pmatrix} \Rightarrow 1=47(-7)+30(11)$

代码

int gcdEx(int a, int b, int *x, int *y) 
 {
     if(b==0)
     {
         *x = 1,*y = 0;
         return a;
     }
     else
     {
         int r = gcdEx(b, a%b, x, y); /* r = GCD(a, b) = GCD(b, a%b) */
         int t = *x;
         *x = *y;
         *y = t - a/b * *y;
         return r;
     }
 }