CodeIQ MAGAZINECodeIQ MAGAZINE

数学の問題をプログラミングで解こう!「ディビジョン・サム」問題解説

2016.12.02 Category:CodeIQ問題解説・リーダーボード Tag: ,

  • 7
  • このエントリーをはてなブックマークに追加

「数学の問題をプログラミングで解こう!」シリーズ。
非常に大きな数の約数の和を求めるという問題でした。
値そのものを計算しないで、約数の和のみをうまく工夫して求めましょう。

みなさんはうまくできましたか?
ぜひ出題者のKawazoeさんによる本記事で解法をチェックして下さい!
by CodeIQ運営事務局

問題のおさらい

自然数 n に対し、(n!)n (つまり n の階乗の n 乗)の約数の和を F(n) と定義します。

例えば、F(2) = 7 です。(2!)2 = 4 で、4 の約数は 1, 2, 4 だからです。
また、F(3) = 600 です。(3!)3 = 216 で、216 には 16 個の約数があり、それらの和は 600 です。
同様に、F(4) = 991111 となります。

標準入力から、自然数 n(1 ≦ n ≦ 103)が与えられるとき、標準出力に F(n) を 1000003 で割った余り を出力するプログラムを書いてください。
たとえば、F(12) を 1000003 で割った余りは 845575 となります。

方針:素因数分解された形で考えよう

(n!)n という非常に大きな数の約数に関する問題です。(n!)n の値を実際に計算してみるという方法は非効率だと多くの方が気づかれたことと思います。

ではどうすればよいでしょうか。(n!)n素因数分解された形で考える ことがポイントです。下記の Wikipedia の記事が参考になります。一般に、整数 k が 、互いに異なる素数 p1, p2, …, pm と、その指数 a1, a2, …, am を使って、p1a1×p2a2× … ×pmam と表されるとき、k の約数の和は次のように表されます。

divisionsum-1

参考:約数:約数の和(Wikipedia)

したがって本問では、(n!)n を素因数分解したとき、上記の p1, p2, …, pma1, a2, …, am がどのような値になるかを考えましょう。もちろん、(n!)n の値を直接計算するのでなく、(n!)n という形に着目して考えましょう。

(n!)n という数は、整数の積として見ると、1 から n までの整数が n 個ずつ集まったものの積です。そこでまずは 1 から n までの整数を、一つずつ素因数分解しましょう。各素数が何回現れるかをそれぞれ求め、n までの総数を計算します。そしてその数を n 倍したものが、(n!)n の素因数の指数と分かります。

以上の方法で、(n!)n を素因数分解したときの p1, p2, …, pma1, a2, …, am の値が求められます。あとは、上の約数の和の公式により、答えを計算できます。

コード例は以下です(Python)。前半では、上記の a1, a2, …, am に対応する配列 a を用意し、素数 p に対する指数を a[p] に保存します。後半では、公式の計算を行っています。公式の各カッコの中は、公比 p の等比数列です。1 から始めて、次々と p をかけながら和をとっていきましょう。以下は Python のコード例です。

n = input()
a = [0] * (n+1)
D = 1000003
for k in range(2, n+1):
    p = 2
    while k != 1:
        while k % p == 0:
            a[p] += n
            k /= p
        p += 1

answer = 1
for p in range(2, n+1):
    t, s = 0, 1
    for j in range(a[p]+1):
        t = (t + s) % D
        s = (s * p) % D
    answer = (answer * t) % D
print answer

このコードで、n ≦ 1000 という条件を 1 秒の制限時間でクリアできます。

さらなる高速化! 等比数列の和の公式を使おう

上のコードでひとまずクリアは可能ですが、実は本問にはさらなる高速化のポイントがあります。注目すべきは後半です。ここでは 1+pp2+…+pa の値を求めるために計 a 回の乗算と加算を行っています。n ≦ 1000 という条件だと、特に小さい pa はけっこう大きな値になります(p=2 なら a は最大で約 106。)

そこで役に立つのが、等比数列の和の公式です。

divisionsum-2

この考えは過去の記事でも登場しました。(参考:数学の問題をプログラミングで解こう!「エース・ナンバー」問題解説の方法2。)「pa+1 乗の計算」と「p-1 での割り算」を行うことがポイントです。それぞれ「べき剰余」と「100003 を法とする p-1 の逆元」を求めることで、上よりもずっと短い計算時間で求められます。

べき剰余と逆元については、下記のサイトを参考にして下さい。

参考:冪剰余(Wikipedia)

参考:整数の合同と1次合同式(外部サイト)

さらに前半で、1 から n の値を素因数分解する部分でも、高速化が可能です。上のコードでは、ある整数 k に対し、割る数を 1 ずつ変えながら割り算を行っていました。ここでは、素数のそれぞれが 1 から n の積に素因数として何個入っているかを調べましょう。たとえば n が 14 の場合、素数 2 が何個入っているかは、図のようにして求めることができます。

divisionsum-3

始めにエラトステネスの篩で素数の一覧を作り、素数のそれぞれについて上記の計算を行いましょう。それを n 倍したものが、(n!)n の各素因数の指数となります。

コード例を以下に示します。はじめのコードでは、n=1000 のとき 0.5 秒程度の実行時間で答えを出せていましたが、今回のコードではそれよりもずっと短い実行時間(実行くんで 0.01 秒程度)で答えを出せます。

# 拡張ユークリッドの互除法
# a*x + b*y = gcd(a,b) となる (x,y) を求める.
def gcd_ext(a, b):
    x, y, lastx, lasty = 0, 1, 1, 0
    while b != 0:
        q = a/b
        a, b = b, a%b
        x, y, lastx, lasty = lastx-q*x, lasty-q*y, x, y
    return (lastx, lasty)

# d を法とする a の逆元 (a*x = 1 (mod d) となる x) を求める.
def inv_mod(a, d):
    (p, q) = gcd_ext(a, d)
    return p % d

# a^n mod d を求める.
def pow_mod(a, n, d):
    if n == 0: return 1
    if n % 2 == 0: return pow_mod((a*a)%d, n/2, d)
    return a * pow_mod(a, n-1, d) % d

n = input()
# エラトステネスの篩を作る
prime = [True] * (n+1)
for i in range(2, n+1):
    for j in range(2*i, n+1, i):
        prime[j] = False

a = [0] * (n+1)
D = 1000003
for p in range(2, n+1):
    if prime[p] == False: continue
    k = p
    while k <= n:
        a[p] += n*(n/k)
        k *= p

answer = 1
for i in range(2, n+1):
    t = (pow_mod(i, a[i]+1, D) - 1) * inv_mod(i-1, D)
    answer = (answer * t) % D
print answer

みんなのコード

本問の挑戦者で、コードを公開して頂いた方のコードを Togetter でまとめました。(随時追加していきます。)こちらもぜひチェックして下さい。

CodeIQ「ディビジョン・サム」 みんなのコード

CodeIQ運営事務局より

kawazoeさん、ありがとうございました!
kawazoeさんの次回の問題にご期待下さい。

  • 7
  • このエントリーをはてなブックマークに追加

■関連記事

【息抜き】カードを上手く並べよう【言語不問】解答と解説... 【息抜き】カードを上手く並べよう【言語不問】 本問題は、表題のテーマで、簡単なプログラムを書くものです。 それでは以下、問題とその解答を見ていきましょう。 問題 あなたは、11から99までの、89枚のカードを持っています。問題では、横幅と高さの整数が与えられます。この横幅と高さで作られるマス...
【コードミステリ】数字に隠されたメッセージ【言語不問】解答と解説... 【コードミステリ】数字に隠されたメッセージ【言語不問】 本問題は、表題のテーマで、簡単なプログラムを書くものです。 喜屋武ちあきさんによるCodeIQ MAGAZINEでのブックレビューに合わせて、『顔貌売人』(文藝春秋)とのコラボ問題として出題されたものです。 それでは以下、問題とその解...
数学の問題をプログラミングで解こう!「ディバイド・アウト」問題解説... 問題のおさらい 自然数 n と素数 p に対し、n の階乗(n!)を p でこれ以上割り切れなくなるまで繰り返し割り、その商をさらに p で割ったときの余りを F(n, p) と定義します。 例えば F(12, 5)=4 です。 12!(=479001600)は 5 で最大 2 回割ることができ...
【息抜き】ファイル名を作ろう【言語不問】解答と解説... 【息抜き】ファイル名を作ろう【言語不問】 本問題は、表題のテーマで、簡単なプログラムを書くものです。 それでは以下、問題とその解答を見ていきましょう。 問題 ファイルをディレクトリ内に作成する際、同じ名前のファイルがあると、末尾に数字を付けるなどして同じ名前にならないようにします。 こうし...
【夏のミステリー】殺人現場のコード 解答と解説... 【夏のミステリー】殺人現場のコード 本問題は、表題のテーマで、簡単なプログラムを書くものです。 それでは以下、問題とその解答を見ていきましょう。 問題 殺人現場にプログラマが倒れていて、途中までプログラムが書かれている。 「続きを書いて欲しい」 これはダイイングメッセージなのか? どう...
数学の問題をプログラミングで解こう!「キャンディ・アンド・チョコレート」問題解説... 問題のおさらい n 個のキャンディをグループに分けます。 グループの最大のキャンディの個数が k 個となるような分け方の数を F(n, k) と定義します。 例えば、F(8, 3)=5 です。このときの分け方を以下に示します。 なお個々のキャンディを区別せずに扱う点に注意してください。 同...

今週のPickUPレポート

新着記事

週間ランキング

CodeIQとは

CodeIQ(コードアイキュー)とは、自分の実力を知りたいITエンジニア向けの、実務スキル評価サービスです。

CodeIQご利用にあたって
関連サイト
codeiq

リクルートグループサイトへ