L 记事本
扩展欧几里得算法
前面的笔记中,ECDH,ECDSA, RSA等整数模素数的操作,都需要用到模逆元
本文介绍 扩展欧几里德算法快速求模逆元
\[e d = 1\pmod n\]裴蜀定理(贝祖等式):
若$a,b$是整数,且$gcd(a,b)=d$，那么对于任意的整数$x,y$,$ax+by$都一定是$d$的倍数，特别地， 一定存在整数$x,y$，使$ax+by=d$成立。
推论:a ,b 互质的充要条件就是 $ax +by = 1$
欧几里得算法求最大公约数
$d$ 是 a,b 的最大公约数(Greatest Common Divisor )
由 $ax +by = kd$ 知道,我们通过选取不同x,y可以逐步让k小,一直到0,那么前一次 肯定就是$k = 1$,
可以让大数减去小的,获得的差取代大数,然后在不停迭代。
例如 a = 78 = 3 * 13 * 2 ,b= 28 = 2 * 7
78 - 28 = 50 , a= 50 , b = 28
50 - 28 = 22 , a = 22 ,b = 28
28 - 22 = 6 , a = 28,b = 6
28 - 6 * 4 = 4 , a = 4 ,b = 6 /// 这里连续减法,重复步骤直接用乘法
6 - 4 = 2 , a = 2 ,b = 4
4 - 2 * 2 = 0 // 上一步的余数 2 就是最大公约数
扩展欧几里得算法
由于a,b 互质.
$ax + by = 1$
这里假定 $a < b$
我们正常欧几里得算法过程
\(k_0 = a , \\ - k_0x_1 + b = k_1 \\ - k_1 x_2 + k_0 = k_2 \\ - k_2 x_3 + k_1 = k_3 \\ - k_3 x_4 + k_2 = k4 \\ \vdots \\ - k_i x_{i+1} + k_{i-1} = k_{i+1} \\ \quad k_i < k_{i-1 } , x_{i+1} = \lfloor{k_{i+1}/k_{i}}\rfloor \\ \vdots \\ - k_n x_{n+1} + k_{n-1} = 1 ,这里可以终止,ab 互质,最大公约数一定是1\) 构造 $ax + by = 1$
思路 1(走不通)
注意到,倒数第1步中有$\dots =1$的等式, \(- k_n x_{n+1} + k_{n-1} = 1\) 将$k_{n-1} = \dots$ 依次代入. 很遗憾,这样会有另外一个常数,会得到 $ax + by + C = 1$ 形式的方程，这个不是我们想要的。
思路2
类似于求斐波那契数列构造
\[\begin{bmatrix} F_n \\ F_{n-1} \end{bmatrix} = A \begin{bmatrix} F_{n-1} \\ F_{n-2} \end{bmatrix}\]由上面的递推公式，
\[\begin{bmatrix} k_{n} \\ k_{n-1} \end{bmatrix}\] \[- k_{n-1} x_{n} + k_{n-2} = k_{n}\] \[\downarrow \\ \begin{bmatrix} k_{n} \\ k_{n-1} \end{bmatrix} = \begin{bmatrix} -x_n & 1 \\ \\ 1 & 0 \end{bmatrix} \begin{bmatrix} k_{n-1} \\ \\ k_{n-2} \end{bmatrix}\]逐步展开,一直到 $b, a$ ,这里 $a = k_0$,$k_n = 1$ (互质,取到最大公约数)
\[\begin{bmatrix} 1 \\ k_n-1 \end{bmatrix} \begin{bmatrix} k_n \\ k_n-1 \end{bmatrix} = \begin{bmatrix} -x_n & 1\\ 1 & 0 \end{bmatrix} \begin{bmatrix} -x_{n-1} & 1 \\ 1 & 0 \end{bmatrix} \dots \begin{bmatrix} -x_{1} & 1 \\ 1 & 0 \end{bmatrix} \begin{bmatrix} a\\ b \end{bmatrix}\]这样,我们迭代欧几里德算法时候,同步算出矩阵的乘积.
到最后,矩阵的 $m_{11} = x$ 是模逆元, $m_{12} = y$ 是方程的根, \(ax + by = 1\) $x 为 a 的模b 逆元$
下面是一个demo代码. 参数内建都是Number类型,大小有限制,作为原理演示。
function gcd(a,b,r0){
var r = b % a ;
if (r == 0) {
return a;
}
return gcd(r,a);
}
function exGcd(a,b){
/**
* 2x2 矩阵,用数组表示
* @param {*} a
* @param {*} b
*/
function matrixMultiply(A,B){
return [A[0] * B[0] + A[1] * B[2],A[0] * B[1] + A[1] * B[3],
A[2] * B[0] + A[3] * B[2] , A[2] * B[1] + A[3] * B[3]];
}
function logN(a,msg){
console.log('' + msg,a.toString(16))
}
function _exGcd(a,b){
/**
* | xn 1 | ... | x1 1| * | a | = | 1 |
* | 1 0 | | 1 0| | b | = | k |
*/
var s = b % a ;
var q = (b - s ) /a ;
if(s == 1){ //终止
return [-q,1n, 1n,0n]
}
else if(s == 0){
/// error 不互质,除非 a == 1
return [0n,0n,0n,0n];
}
const matrix0 = [-q,1n ,1n,0n];
const matrix1 = _exGcd(s,a);
return matrixMultiply(matrix1,matrix0);
}
if(a > b) {
a = a % b ;
}
/// 1 的逆元是本身
if (a == 1) {
return 1n;
}
// 0 无逆元
if (a == 0) {
return 0;
}
var m = _exGcd(a,b);
var r = (m[0] % b + b) % b
return r
}
// var b = 10;
// for(var a = 1 ;a < b ; ++ a ){
// var c = exGcd(a,b);
// if (c == 0 || c == 1) {
// // console.log('-----' ,a)
// }
// else if(a * c %b == 1){
// console.log(c , a , ` ${c} * ${a} % ${b} = ${a * c %b }`);
// }
// }
module.exports = exGcd
二进制扩展欧几里得算法
这里还是用到了很多乘法，在数字很大情况下，效率低下。