Skip to content

Instantly share code, notes, and snippets.

@kuuso
Created March 10, 2017 15:31
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save kuuso/00235eacbe29f5e89a738937a769f43f to your computer and use it in GitHub Desktop.
Save kuuso/00235eacbe29f5e89a738937a769f43f to your computer and use it in GitHub Desktop.
using System;
using System.Collections;
using System.Collections.Generic;
using System.Linq;
using System.Text;
class TEST{
static void Main(){
Sol mySol =new Sol();
mySol.Solve();
}
}
class Sol{
public void Solve(){
// G(n) := F(n!) とする
#if DEBUG
// sample
Check(G(10), 829440);
Check(G(20), 961998);
Check(G(1), 1);
Check(G(2), 1);
Check(G(3), 2); // F(6) = 2 ({1,5})
Check(G(4), 8); // F(24) = 8 ({1,5,7,11,13,17,19,23})
#endif
int N = ri();
Console.WriteLine(G(N));
}
long G(int N){
// F(n) = n / Π p_i * Π (p_i - 1) where {p_i} prime, p_i < n for all i
// 素数で割る部分は変わらず,先頭を n! に置き換えればいいだけでは.
long mod = 1000003;
bool[] isPrime = new bool[N+1];
for(int i=2;i<N+1;i++) isPrime[i] = true;
long ret = 1;
for(int i=2;i<=N;i++){
ret *= i; ret %= mod;
if(!isPrime[i]) continue;
ret *= (i-1); ret %= mod;
ret *= ModInv(i, mod); ret %= mod;
for(int j=i+i;j<=N;j+=i) isPrime[j] = false;
}
return ret;
}
public static long ModInv(long x, long mod){
//拡張ユークリッドの互除法でmod逆元を求める
//ax + mod * y = 1 で mod取れば ax = 1 になるから xがaの逆元になる
//逆元が存在しない→0を返す。
long a = 0, b = 0, c = 0;
if(x == 0) return 0;
ExtGCD(x, mod, ref a, ref b, ref c);
if(c != 1) return 0;
return (a + mod) % mod;
}
public static void ExtGCD(long x,long y,ref long a,ref long b,ref long c){
long r0 = x; long r1 = y;
long a0 = 1; long a1 = 0;
long b0 = 0; long b1 = 1;
long q1, r2, a2, b2;
while(r1 > 0){
q1 = r0 / r1;
r2 = r0 % r1;
a2 = a0 - q1 * a1;
b2 = b0 - q1 * b1;
r0 = r1; r1 = r2;
a0 = a1; a1 = a2;
b0 = b1; b1 = b2;
}
c = r0;
a = a0;
b = b0;
}
static bool Check(long input, long ans){
Console.Write(input==ans?"OK":"NG");
if(input!=ans){
Console.Write(": input={0},ans={1}",input,ans);
}
Console.WriteLine();
return input==ans;
}
static int ri(){ return int.Parse(Console.ReadLine()); }
}
totient関数φ(n)をつかいます.問題のF(n) は φ(n) そのものです.
totient関数φ(n)を求める方法として,nの素因数分解 n = Π p_i に関して
φ(n) = n / ( Π (p_i - 1)) * Π p_iという公式があります.
n!の素因数は nの素因数と一致するので,結局 n!を求める途中で見つかった素数pに関して
* (p-1) / p を実施してやればよいです.
コードではエラトステネスの篩を同時に行いながら素数を見つけては* (p-1) / p をmod 1000003 演算で実施しています.
(mod逆元を使わなくても,階乗計算時にpを掛けるのを省略するだけで良かったな...)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment