はじめに
- 初級編では、nが10000程度のケースまでしか対応できなかったので、10^7程度まで対応できるような実装を紹介します。
- ただし、ここから先は、modについて多少の知識が必要となります。具体的には、(mod M)における逆元を使用します。
- だから、この投稿は、ほぼ、(mod M)の逆元の話になるよ。
(mod M)における逆元について
- 例えば、整数aに対する「aの逆元」を考えます。$a^{-1}$と書き、ax ≡ 1 (mod M) を満たすxの解です。
- 詳細は、過去の投稿を参考して欲しいけど、(mod 10)における逆元の具体的なイメージを書けば、
- $1^{-1}$ ≡ 1 、言うまでもなく、全てのmに対して、1×1 ≡ 1 (mod m)となる。
- $2^{-1}$ は存在しない。法10と2は素でないため、2x ≡ 1 を満たす解xは存在しない。感覚的にも、2にどんな数を乗じても、下一桁が1になることがないと分かるよね。
- $3^{-1}$ ≡ 7 、これは 3×7=21 より分かる。
- $4^{-1},5^{-1},6^{-1}$ は存在しない。
- $7^{-1}$ ≡ 3 、これは 7×3=21 より分かる。
- $8^{-1}$ は存在しない。
- $9^{-1}$ ≡ 9 、これは 9×9=81 より分かる。
- $10^{-1}$ は存在しない。
- $11^{-1}$ ≡ 1 、これは 11×1=11 より分かる。
どうして逆元で解決できるの？
- 初級編で書いた脊髄反射で書いたような下記コードでは、②部分で割り算が登場したので、mod化出来なかったね。
func comb(_ a:Int,_ b:Int)->Int{
var c = a+b
var a = a
var b = b
var ans = 1
while c>1 { ans *= c ; c-=1 } // ①
while a>1 { ans /= a ; a-=1 } // ②
while b>1 { ans /= b ; b-=1 }
return ans
}
- だけど、例えば(mod m)におけるaの逆数をinvMOD(a,m)とする関数invMODを用意できれば、②は
while a>1 { ans *= invMOD(a,m) ; a-=1 }と書き換えることが出来るため、mod化が楽勝で出来るようになるね。 - そうすれば、高速化も多少は期待できるかな？まぁ、invMOD関数自体の計算量が重いと駄目だけどね。
逆元を算出する方法
- まず、ベズーの等式を紹介します。
ベズーの等式
a と b を 0 でない整数とし、gcd(a,b) をそれらの最大公約数とする。このとき、以下の等式を成立させる整数 x と y が存在する。
ax + by = gcd(a,b)
※gcd:最大公約数(greatest common divisor)
- つぎに、次の等式を紹介します。
ベズーの等式「っぽい見た目の」等式
\displaylines{
n･n^{-1} + mk ≡ 1 \pmod m
}
※ nとmは素 ( ∵素でないと、nの逆元が存在しない)
※ $n･n^{-1}≡1$より、上記等式は任意のkについて成立する。
- ここで、互いに素なnとmについて、ベズーの等式を当てはめると、gcd(n,m)=1より以下のことが言える。
-
nx+my = 1を成立させるx,yが存在する。
-
- よって、上記により求めた(x,y)を$(n^{-1},k)$と考えることが出来る。
- このように、$n^{-1}(=x)$を求める事ができます。
- なお、(a,b)より(x,y)を求めるコードは以下のとおりです。コードの解説は、こちら。
func exGCD(_ a:Int,_ b:Int)->(Int,Int){
if (b == 0) {return (1,0)}
let temp = exGCD(b,a % b)
let q = a / b
return (temp.1 , temp.0 - q * temp.1)
}
let a = 32
let b = 20
let ans = exGCD(a,b)
print( ans.0, ans.1, a*ans.0 + b*ans.1 ) // 2 -3 4
- このexGCDを利用して、逆元を求めるコードは、以下のとおり。
func invMOD(_ n:Int,_ mod:Int)->Int{
var inv = exGCD(n,mod).0
if inv < 0 {inv += mod}
return inv
}
脊髄反射コードに逆元を導入してmod化する
- もう、大した事は残ってないよ。
func comb(_ a:Int,_ b:Int)->Int{
let mod = 1000000007
var c = a+b
var a = a
var b = b
var ans = 1
while c>1 { ans = ans * c % mod ; c-=1 } // ①
while a>1 { ans = ans * invMOD(a,mod) % mod ; a-=1 } // ②
while b>1 { ans = ans * invMOD(b,mod) % mod ; b-=1 }
return ans
}
func invMOD(_ n:Int,_ mod:Int)->Int{
var inv = exGCD(n,mod).0
if inv < 0 {inv += mod}
return inv
}
func exGCD(_ a:Int,_ b:Int)->(Int,Int){
if (b == 0) {return (1,0)}
let temp = exGCD(b,a % b)
let q = a / b
return (temp.1 , temp.0 - q * temp.1)
}
print(comb(10_000_000,10_000_000)) // 486682686
- 上記のとおり、$10^7$のオーダーの計算は出来ました。ちなみに$2×10^7$にしたらpaizaではtimeoutでした。
おわりに
- ほとんど、過去の投稿を再構成しただけだったね。
- でもこれで、桁が大きい場合の計算が出来るようになったよ！