MOD 逆数 - 幸せに

$n$ を整数とすると， $n$ と互いに素な任意の整数 $a$ に対し， $ax \equiv 1\ ({\rm mod}\ n)$ を満たす整数 $x\ (0 \leq x < n)$ がただ一つ存在します．この $x$ を， $n$ を法としたときの逆数といい，競技プログラミングで頻繁に役立ちます．今回はこれを計算するコードを C で書いていきます．

数学的な解き方

$ax \equiv 1\ ({\rm mod}\ n)$ ならば，ある整数 $x'$ が存在して $ax + nx' = 1$ です．すると，これは $nx' \equiv 1\ ({\rm mod}\ a)$ と書き換えられます．ここで $n$ を $a$ で割った余りを $r$ とすると， $rx' \equiv 1\ ({\rm mod}\ a)$ となります．これを繰り返すことで $x$ を求めることができます．具体的には， $n = 529, \, a = 100$ とすると

$100x \equiv 1\ ({\rm mod}\ 529)$ ある $x'$ が存在して $100x + 529x' = 1$ よって $529x' \equiv 1\ ({\rm mod}\ 100)$
$29x' \equiv 1\ ({\rm mod}\ 100)$ ある $x''$ が存在して $29x' + 100x'' = 1$ よって $100x'' \equiv 1\ ({\rm mod}\ 29)$
$13x'' \equiv 1\ ({\rm mod}\ 29)$
$13x'' + 29x''' = 1$
$29x''' \equiv 1\ ({\rm mod}\ 13)$
$3x''' \equiv 1\ ({\rm mod}\ 13)$
$3x''' + 13x'''' = 1$
$13x'''' \equiv 1\ ({\rm mod}\ 3)$
$x'''' \equiv 1\ ({\rm mod}\ 3)$
これより，最後の $x''''$ の値を $1$ とすれば，逆に辿って $x$ の値を求めることができます．逆に辿るやり方は，
$\begin{align*}3x''' &= 1-13x'''' \\ &= 1 - 3\cdot 4x'''' - x'''' \\ &= - 3\cdot 4x''''\end{align*}$ より
$x''' = -4 x'''' = -4$
$\begin{align*}13x'' &= 1 - 29x''' \\ &= 1 - 13\cdot 2x''' - 3x''' \\ &= 1 - 13\cdot 2x''' - (1 - 13 x'''') \\ &= -13\cdot 2x''' + 13 x''''\end{align*}$ より
$x'' = -2x''' + x'''' = 9$
$\begin{align*}29x' &= 1 - 100x''\\ &= 1 - 29\cdot3x''-13x''\\ &= 1 - 29\cdot3x''-(1-29x''')\\ &= -29\cdot3x'' + 29 x'''\end{align*}$ より
$x' = -3x'' + x''' = -31$
$\begin{align*}100x &= 1 - 529x' \\ &= 1 - 100\cdot5x' - 29x' \\ &= 1 - 100\cdot 5x' - (1 - 100 x'') \\ &= - 100 \cdot 5x' + 100 x''\end{align*}$ より
$x = -5x' + x'' = 164$ という具合です．

一般化

上の例より， $a_0 = n$ ， $a_1 = a$ とし，

$\begin{align*} a_0 \div a_1 &= b_0 \cdots a_2 \\ a_1 \div a_2 &= b_1 \cdots a_3 \\ & \vdots \\ a_m \div a_{m + 1} &= b_m \cdots 1 \end{align*}$ とすると， $\begin{align*} x_{m + 1} &= 1 \\ x_m &= - b_m x_{m + 1} \\ x_{m - 1} &= - b_{m - 1} x_m + x_{m + 1} \\ x_{m - 2} &= - b_{m - 2} x_{m - 1} + x_m \\ & \vdots \\ x_0 &= -b_0 x_1 + x_2 \end{align*}$ という手続きで $x = x_0$ の値が求められることが分かりました．

行列

上の手続きを行列を用いて表すと

$\begin{pmatrix} x_k \\ x_{k + 1} \end{pmatrix} = \begin{pmatrix} -b_k & 1 \\ 1 & 0 \end{pmatrix} \begin{pmatrix} x_{k + 1} \\ x_{k + 2} \end{pmatrix}$ (ただし $k$ は $0$ 以上 $m$ 以下の整数，また $\begin{pmatrix}x_{m + 1} \\ x_{m + 2}\end{pmatrix} = \begin{pmatrix}1 \\ 0 \end{pmatrix}$ とする)となります．つまり $\begin{pmatrix} x_0 \\ x_1 \end{pmatrix} = \begin{pmatrix} -b_0 & 1 \\ 1 & 0 \end{pmatrix} \begin{pmatrix} -b_1 & 1 \\ 1 & 0 \end{pmatrix} \cdots \begin{pmatrix} -b_m & 1 \\ 1 & 0 \end{pmatrix} \begin{pmatrix} 1 \\ 0 \end{pmatrix}$ なので，行列計算によって $x_0$ を求めることができます．

コード

この行列計算をそのままコードで書くと次のようになります．

int inverse(int mod, int val){
	int x[2] = {mod, val};
	int a[2][2][2] = {{{1, 0}, {0, 1}}};
	int i;
	for(i = 0; x[!i]; i ^= 1){
		a[!i][0][0] = a[i][0][1] - x[i] / x[!i] * a[i][0][0];
		a[!i][0][1] = a[i][0][0];
		a[!i][1][0] = a[i][1][1] - x[i] / x[!i] * a[i][1][0];
		a[!i][1][1] = a[i][1][0];

		x[i] = x[i] % x[!i];
	}
	return a[!i][0][0];
}

実際には行列の第 1 行と第 2 行の計算は独立しているので次のように書けます．

int inverse(int mod, int val){
	int x[2] = {mod, val};
	int a[2][2] = {{1, 0}};
	int i;
	for(i = 0; x[!i]; i ^= 1){
		a[!i][0] = a[i][1] - x[i] / x[!i] * a[i][0];
		a[!i][1] = a[i][0];

		x[i] = x[i] % x[!i];
	}
	return a[!i][0];
}

さらに第 1 列と第 2 列を配列内で毎回ひっくり返すようにすると次のようになります．

int inverse(int mod, int val){
	int x[2] = {mod, val};
	int a[2] = {1, 0};
	int i;
	for(i = 0; x[!i]; i ^= 1){
		a[!i] -= x[i] / x[!i] * a[i];
		x[i] = x[i] % x[!i];
	}
	return a[!i];
}

ただしこれだと答えが正になったり負になったりするので，常に $0 \leq x < n$ の範囲で値を得たいときは次のようにします．

int inverse(int mod, int val){
	int x[2] = {mod, val};
	int a[2] = {1, 0};
	int i;
	for(i = 0; x[!i]; i ^= 1){
		a[!i] -= x[i] / x[!i] * a[i];
		x[i] = x[i] % x[!i];
	}
	if(!i) a[!i] += mod;
	return a[!i];
}

この記事があなたの役に立てば幸いです．