「プライム・ペア」問題

「プライム・ペア」問題

締め切りが 2017/03/09 10:00 AM なので，その 1 分後に投稿されるように予約

設問

自然数 k に対し、1 から k までの自然数のうち k と互いに素なものの個数を F(k) と定義します。
（F(k) はオイラーのΦ関数とも呼ばれています。参考：オイラーのφ関数（Wikipedia））
例えば F(12)＝4 です。1 から 12 のうち 12 と互いに素なのは 1, 5, 7, 11 の 4 つです。
標準入力から、自然数 n（1 ≦ n ≦ 105）が与えられます。
標準出力に F(n！) を 1000003 で割った余りを出力するプログラムを書いてください。
例えば n＝10 のとき、F(10！)＝F(3628800)＝829440 です。
同様に、F(20！) を 1000003 で割った値は 961998 です。

===================================

F = function(n) {
    mx = n # n 以下の素数リストを prime に作る
    tbl = 1:mx
    tbl[1] = 0
    for (i in 2:floor(sqrt(mx))) {
        if (tbl[i]) {
            mx2 = mx %/% i
            tbl[2:mx2*i] = 0
        }
    }
    prime = tbl[tbl > 0]
# F(n!) = n! * Π(1-1/prime[i]) = 1*2*3*…*n × (prime[1]-1)/prime[1] × (prime[2-1)/prime[2] …
# F(20!) = 1*2*3*4*5*6*7*8*9*10*11*12*13*14*15*16*17*18*19*20 * (1-1/2)*(1-1/3)*(1-1/5)*(1-1/7)*(1-1/11)*(1-1/13)*(1-1/17)*(1-1/19)
#        = 1*2*3*4*5*6*7*8*9*10*11*12*13*14*15*16*17*18*19*20 * 1/2 * 2/3 * 4/5 * 6/7 * 10/11 * 12/13 * 16/17 * 18/19
#        = 1*4*6*8*9*10*12*14*15*16*18*20 * 1*2*4*6*10*12*16*18
#        = 416084687585280000
# 416084687585280000 %% 1000003 = 961998
# i = 1 ~ n までの積。ただし，i が素数の場合は「i-1」
    x = 1:n
    x[prime] = x[prime] - 1
    ans = 1
    for (a in x) { # 要素の掛算（ただし，毎回 1000003 の剰余を結果とする）
        ans = (ans*a) %% 1000003
    }
    cat(ans)
}
#F(scan(file("stdin", "r")))

F(10) # 829440
F(20) # 961998
F(30) # 845477
F(100) # 145380
F(431) # 945906
F(2017) # 272985
F(81645) # 603326
F(99100) # 799614