hdu5407 CRB and Candies

令f(n) = lcm((n, 0), (n, 1),..., (n, n))，

g(n) = lcm(1, 2,.., n)。

则有结论：f(n) = g(n + 1) / (n + 1)（证明过程参阅：http://arxiv.org/pdf/0906.2295v2.pdf）。

待求的是ans = f(n) % mod

= (g(n + 1) / (n + 1)) % mod

= (g(n + 1) % mod * Inverse(n + 1) % mod) % mod。

其中Inverse(n)表示n的关于mod的乘法逆元，即满足Inverse(n) * n % mod = 1。

由费尔马小定理，a^mod-1%mod = 1， 于是 a * a^{mod -2} % mod = 1。

则Inverse(n) = a^mod-2 %mod(mod 为素数），逆元唯一。

这样用快速幂很快计算一个数的乘法逆元，下面需要计算g(n)。

有以下结论g(n) = 

g(n - 1) * p (如果存在质数p和整数k，满足n = p^k)（1），

g(n - 1) (否则）。

我们需要预处理出一个对于每个整数其可能对应的p的数组以快速计算g(n)。

所有的素数存在唯一的p等于其本身和k=1满足式（1）。

对于任一整数可以写成素数幂的乘积的形式，因此可以用当前素数往后更新。

http://acm.hdu.edu.cn/showproblem.php?pid=5407

#include <cstdio>
#include <cstring>
#include <algorithm>

using namespace std;
typedef __int64 LL;

const int maxn = 1e6 + 10;
const int mod = 1e9 + 7;

int n;
LL f[maxn], g[maxn];

LL power(LL a, LL p){
    LL ans = 1;
    a %= mod;
    while(p){
        if(p & 1){
            ans = (ans * a) % mod;
        }
        p >>= 1;
        a = (a * a) % mod;
    }
    return ans;
}

LL Inverse(LL num){
    return power(num, mod - 2);
}

void solve(){
    LL ans = f[n + 1] * Inverse(n + 1) % mod;
    printf("%I64d
", ans);
}

bool check(int num){
    int d = g[num];
    while(num % d == 0 && num > 1) num /= d;
    return num == 1;
}

void init(){
    for(int i = 1; i < maxn; i++) g[i] = i;
    for(int i = 2; i < maxn; i++){
        if(g[i] == i){
            for(int j = i + i; j < maxn; j += i){
                g[j] = i;
            }
        }
    }
    f[0] = 1;
    for(int i = 1; i < maxn; i++){
        if(check(i)) f[i] = f[i - 1] * g[i] % mod;
        else f[i] = f[i - 1];
    }
}

int main(){
    init();
    int T;
    scanf("%d", &T);
    while(T--){
        scanf("%d", &n);
        solve();
    }
    return 0;
}