LoginSignup
0
1

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

@kei_tnk

二項係数をmod表示する（上級編）【競技プログラミング】

Posted at

はじめに

  • 中級編では、$10^7$のオーダーで二項係数を計算できるようにしたよ。
  • でも、速度的に問題が･･･
  • よって、今回は高速化を目指すよ。
  • その前に、modに関する新知識「フェルマーの小定理」も紹介しちゃいます。

フェルマーの小定理

pが素数であるとき、下記が成立する。

a^p ≡ a \pmod p

さらに、pとaが互いに素(gcd(p,a)=1)であるとき、下記が成立する

a^{p-1} ≡ 1 \pmod p
  • 上記の証明は、結構分かり易いから、一読しておくと良いよ。
  • で、既に気がついた人がいるかもしれないけど、これが分かるね。

pが素数であり、pとaが互いに素であるとき、

a･a^{p-2} ≡ 1 \pmod p

であるから、$a^{-1}$は$a^{p-2}$となる

  • 正直言って、「中級編より簡単じゃ〜ん、早く言ってよ〜〜！！」って感じだよね。
  • でもさ、気付いてる？コンテストで出てくるmodって、(mod 1000000007)みたいな感じよ。
  • 累乗計算も馬鹿正直に実行すると、計算量がかかっちゃうから、繰り返し二乗法を使うよ。
func powMOD(_ a:Int,_ n:Int,_ mod:Int)->Int {
    var ans = 1
    var x = a
    var i = n
    while i > 0 {
        if i % 2 == 1 {ans *= x}
        x *= x
        i >>= 1
        
        if ans > mod {ans = ans % mod}
        if x > mod {x = x % mod}
    }
    return ans
}

let mod = 998244353

print(powMOD(11,3,mod)) // 1331
print(powMOD(3,33,mod)) // 478528062
print(powMOD(7,127,mod)) // 561762571

フェルマーの小定理をつかって逆元を求めるコード

  • じゃあ、やってみよ〜
func invMOD(_ n:Int,_ mod:Int)->Int{
    powMOD(n,mod-2,mod)
}
  • ざっと、こんなもんだね！
  • 計算量を軽く検証してみよう！
let mod = 1000000007
var cnt = 0
for i in 1000000000...1000000006 {
    print(i,invMOD(i,mod),cnt)
    cnt = 0
}

func powMOD(_ a:Int,_ n:Int,_ mod:Int)->Int {
    var ans = 1
    var x = a
    var i = n
    while i > 0 {
        cnt += 1 // 計算量チェック用に追加！！！
        if i % 2 == 1 {ans *= x}
        x *= x
        i >>= 1
        
        if ans > mod {ans = ans % mod}
        if x > mod {x = x % mod}
    }
    return ans
}

func invMOD(_ n:Int,_ mod:Int)->Int{
    powMOD(n,mod-2,mod)
}
  • while関数にcnt変数を追加したよ。そして、mod=1000000007より、少し小さい1000000000...1000000006をnとして、逆元を計算してみたよ。その結果は以下のとおりだから、30回転しているね。
1000000000 857142863 30
1000000001 833333339 30
1000000002 600000004 30
1000000003 750000005 30
1000000004 666666671 30
1000000005 500000003 30
1000000006 1000000006 30

中級編の速度を検証してみる

  • フェルマーの小定理を使ったコードと同じようにcnt変数を追記して検証したよ。
let mod = 1000000007
var cnt = 0
for i in 1000000000...1000000006 {
    print(i,invMOD(i,mod),cnt)
    cnt = 0
}

for i in 100000000...100000006 {// 一桁落とした
    print(i,invMOD(i,mod),cnt)
    cnt = 0
}

func exGCD(_ a:Int,_ b:Int)->(Int,Int){
    cnt += 1 // 検証用に追加したよ！！！
    if (b == 0) {return (1,0)}
    let temp = exGCD(b,a % b)
    let q = a / b
    return (temp.1 , temp.0 - q * temp.1)
}

func invMOD(_ n:Int,_ mod:Int)->Int{
    var inv = exGCD(n,mod).0
    if inv < 0 {inv += mod}
    return inv    
}
  • 結果はコチラ！
1000000000 857142863 6
1000000001 833333339 6
1000000002 600000004 6
1000000003 750000005 6
1000000004 666666671 6
1000000005 500000003 5
1000000006 1000000006 4
100000000 571428574 6
100000001 333333339 7
100000002 538461543 8
100000003 391304351 9
100000004 393939397 9
100000005 953488379 8
100000006 566037740 10
  • 回転数は4〜10だから、フェルマーの小定理より回転数は少ないね。ただし、nの桁が大きい方が、回転数が大きいと限らないのも面白いね。

閑話休題

  • 『ちょっと、ちょっと！フェルマーの小定理を使ったのに回転数増えちゃってんじゃん！高速化はどうした！』とお怒りの皆さん。
  • すみません。高速化はこれからです。
  • フェルマーの小定理は、上級編だから一応、手法の取りこぼしのないように紹介しただけでした。
  • でも、高速化って、逆元を求めることを高速化するわけじゃなくて、「二項係数の算出を繰り返す算出する」ことを目指すわけです。
    • もっと具体的に言えば、逆元計算だけでなく、階乗計算も繰り返し計算されることを前提に高速化します。

逆元をまとめて計算する

  • 二項係数を算出するcomb関数が呼び出される前に、逆元計算と階乗計算を予め終えておきたい。
  • メモ化再帰関数とかでなく、最初から、用意しちゃう方法です。
  • だって、二項係数の計算で必須な階乗計算で毎回k!(計算量O(k))を計算するのは馬鹿らしいじゃないですか〜。
    • だから、同じオーダーの計算量O(n)をつかって、1!、2!、3!、･･･、k!、(k+1)! (=k!*(k+1))、･･･、n!を一度に全部計算しちゃっておこうと言うことです。
  • では、逆元も、同じように一度にまとめて計算して、計算量を少なく出来るのだろうかな？....できるんです！

逆元配列inv[i]をボトムアップで求める方法

  1. まず、配列invを用意します。
  2. inv[1] = 1
  3. 法をmとして、i=2...nで、inv[i] = m - (m / i) * inv[m % i] % mを計算します。

※ mを素数とし、n<mとする。よって、i(in 2...n)は全てmと素となる。

  • え？「なんか、気持ち悪い式が出てきた！」って思った？本当に気持ち悪いよね。
  • 気持ち悪いから、解説します。
    スクリーンショット 0006-11-30 22.36.49.jpg
  • 上記より、$i^{-1}$ ≡ -q * inv[m mod i]になります。だけど、inv配列の値は、0..<mの範囲に収めたいから、頭にmを加え、お尻に% mを追加して、以下のようにしたってこと！
    • inv[i] = m - (m / i) * inv[m mod i] % m
  • ご納得いただけただろうか？
  • 結局、「まとめて計算」方式だと、フェルマーの小定理を使った算出方法も、中級編の算出方法も使わなかったね。

高速に二項係数を繰り返し計算する

  • 準備は万端、あとは結果をご覧じろ!
class Comb {
    let mod:Int
    var fac:[Int] // 階乗計算結果
    var inv:[Int] // 逆数計算結果
    var finv:[Int] // 逆数の階乗計算結果
    
    init(mod:Int,max:Int){
        self.mod = mod
        let MAX = max + 1
        
        fac = [Int](repeating: 1, count: MAX)
        inv = [Int](repeating: 1, count: MAX)
        finv = [Int](repeating: 1, count: MAX)       
        
        for i in 2..<MAX {
            fac[i] = fac[i-1] * i % mod
            inv[i] = mod - (mod / i) * inv[mod % i] % mod
            finv[i] = finv[i-1] * inv[i] % mod
        }
        
    }
    
    func comb(_ n: Int, _ k: Int) -> Int {
        fac[n] * finv[k] % mod * finv[n - k] % mod
    }

}

let comb:(Int,Int)->Int = Comb(mod:1000000007,max:20000000).comb

print(comb(12345678,1234567)) // 855321431
  • おお！良い感じにでけた！
  • maxは事前準備しておく階乗計算や逆元計算の数です。$2×10^7$はイケたけど、$3×10^7$だと、paizaではエラーになったよ。
  • あぁ、そう言えば、逆元の階乗計算結果を格納するfinvも導入したよ。

おわりに

おまけ

  • 「逆元をまとめて計算する」の項目で、
    • ＞「まとめて計算」方式だと、フェルマーの小定理を使った算出方法も、中級編の算出方法も使わなかったね。
  • とか書いちゃったけど、実は、一度、フェルマーの小定理等で「n!」の逆元を求めた後、(n-1)!の逆元、(n-2)!の逆元、のように遡って求めていく方法がある。
  • 結論から言えば、iの階乗i!の逆元をfinv[i]に格納するよう事とした場合、以下のように書ける。
    • finv[i] = finv[i+1] * (i+1) % m
    • これは逆向き漸化式って言えば良いのだろうか？
  • なんで、コレが成り立つかを考えれば、こんな感じ。
    • $i!^{-1} × i!　≡　(i+1)!^{-1} × (i+1)!　≡　1\pmod m$
    • このときmが素数で、i<mなら、i!とmは互いに素なので、i!で割っても合同が保たれて、
    • $i!^{-1} ≡　(i+1)!^{-1} × (i+1) \pmod m$
    • お分かりいただけただろうか？
  • だから、Combの一部を下記のように書き換えることが可能
/// 変更前
        for i in 2..<MAX {
            fac[i] = fac[i-1] * i % mod
            inv[i] = mod - (mod / i) * inv[mod % i] % mod
            finv[i] = finv[i-1] * inv[i] % mod
        }

/// 変更後
        for i in 2..<MAX {
            fac[i] = fac[i-1] * i % mod
        }
        finv[MAX - 1] = invMOD(fac[MAX -1])
        for i in (2...MAX-2).reversed() {
            finv[i] = finv[i+1] * (i+1) % mod
        }
  • もちろん、逆元を求める関数invMODは別途用意する必要があるけどね。
    • inv[i] = mod - (mod / i) * inv[mod % i] % modを使うのが気持ち悪いよ〜、とか
    • せっかく、フェルマーの小定理でinvMOD関数を作ったのにもったいないよ〜、とか
    • そんな風に思った人は、こっちで解いても良いかもね。
0
1
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
0
1

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?