danihao123's Blog

有个很好康的结论：

\[\gcd(x^a - 1, x^b - 1) = x^{\gcd(a, b)} - 1\]

然后尝试去用常见套路（枚举gcd）化简柿子，得到：

\[\sum_{k = 1}^n (x^k - 1)\sum_{1\leq a,b\leq n} [\gcd(a, b) = k]\]

推到这里了，看到那个\(\sum_{1\leq a,b\leq n} [\gcd(a, b) = k]\)可能很多同学准备直接上莫比乌斯函数了……其实并不需要，其实那个柿子就是：

\[2\sum_{i = 1}^{\lfloor \tfrac{n}{k}\rfloor} \varphi (i) - 1\]

这个意义还是很显然的……更好的一点是这个柿子可以预处理出来，然后整除分块直接搞一波即可。

代码：

#include <cstdio>
typedef long long ll;
const ll ha = 1000000007LL;
const int maxn = 1000005;
const int N = 1000000;
ll phi[maxn], phi_S[maxn];
void sieve() {
  static bool vis[maxn];
  static int prm[maxn]; int cnt = 0;
  phi[1] = 1LL;
  for(int i = 2; i <= N; i ++) {
    if(!vis[i]) {
      prm[cnt ++] = i;
      phi[i] = i - 1;
    }
    int v;
    for(int j = 0; j < cnt && (v = i * prm[j]) <= N; j ++) {
      vis[v] = true;
      if(i % prm[j] == 0) {
        phi[v] = (phi[i] * (ll(prm[j]))) % ha;
        break;
      } else {
        phi[v] = (phi[i] * phi[prm[j]]) % ha;
      }
    }
  }
  for(int i = 1; i <= N; i ++) {
    phi_S[i] = (phi_S[i - 1] + phi[i]) % ha;
  }
  for(int i = 1; i <= N; i ++) {
    phi_S[i] = ((2LL * phi_S[i]) % ha - 1LL + ha) % ha;
  }
}

ll pow_mod(ll a, ll b) {
  ll ans = 1LL, res = a;
  while(b) {
    if(1LL & b) ans = (ans * res) % ha;
    res = (res * res) % ha;
    b >>= 1;
  }
  return ans;
}
ll inv(ll x) {
  return pow_mod(x, ha - 2LL);
}

ll pre_sum(ll x, ll r) {
  ll p = pow_mod(x, r);
  p = (p - 1LL + ha) % ha;
  p = (p * inv(x - 1LL)) % ha;
  return p;
}
ll seg_sum(ll x, int a, int b) {
  if(x == 1LL) return 0;
  ll len = b - a + 1;
  ll ret = (pre_sum(x, b + 1) - pre_sum(x, a) + ha) % ha;
  ret = (ret - len + ha) % ha;
  return ret;
}
ll calc(ll x, int n) {
  ll ans = 0;
  for(int i = 1; i <= n;) {
    int v = n / i;
    int nx = n / v;
    ll delta = (seg_sum(x, i, nx) * phi_S[v]) % ha;
    ans = (ans + delta) % ha;
    i = nx + 1;
  }
  return ans;
}

int main() {
  sieve();
  int T; scanf("%d", &T);
  while(T --) {
    int x, n; scanf("%d%d", &x, &n);
    printf("%lld\n", calc(x, n));
  }
  return 0;
}

这个题啊……亦可赛艇！

答案是[tex]\varphi(m!)*n!/m![/tex]。原因很简单，把[tex][1,n!][/tex]分成长度为[tex]m![/tex]的若干段，除去第一段外每一段中与[tex]m![/tex]互质的数[tex]k[/tex]肯定满足[tex](k\bmod m!,m!)=1[/tex]（否则，[tex]k[/tex]和[tex]m![/tex]就会有大于一的公因子了)。所以说每一段内与[tex]m![/tex]互质的数都有[tex]\varphi(m!)[/tex]个。

麻烦好像就在于求一个阶乘的欧拉函数。考虑一个新乘上的数能给答案带来的贡献——如果这个数是合数，它的所有因子在前面都有了，只能把他自己贡献出来；如果这个数是质数（假设为[tex]p[/tex]），出了贡献自己以外还会贡献一个[tex](1-1/p)[/tex]，最后乘起来就是贡献了[tex]p-1[/tex]。筛一遍素数再递推一下就好辣~

并且……[tex]n-m[/tex]可能非常大，所以说除去[tex]m![/tex]那块要用逆元做。

（顺便说下阶乘也要递推）

代码：

/**************************************************************
    Problem: 2186
    User: danihao123
    Language: C++
    Result: Accepted
    Time:9408 ms
    Memory:166836 kb
****************************************************************/
 
#include <cstdio>
#include <cmath>
typedef unsigned long long ull;
const int maxn=10000000;
ull R;
bool vis[maxn+5];
inline void sievePrime(){
    register int i,j,m=sqrt(maxn+0.5);
    for(i=2;i<=m;i++)
        if(!vis[i])
            for(j=i*i;j<=maxn;j+=i)
                vis[j]=true;
}
ull fac[maxn+5];
inline void sieveFac(){
    register int i;
    fac[0]=1%R;
    for(i=1;i<=maxn;i++)
        fac[i]=(fac[i-1]*(i%R))%R;
}
ull phifac[maxn+5];
inline void sievePhiFac(){
    register int i;
    phifac[1]=1%R;
    for(i=2;i<=maxn;i++){
        if(vis[i])
            phifac[i]=(phifac[i-1]*(i%R))%R;
        else
            phifac[i]=(phifac[i-1]*((i%R-1%R+R)%R))%R;
    }
}
void exgcd(ull a,ull b,ull& d,ull& x,ull& y){
    if(!b){
        d=a;
        x=1;
        y=0;
    }else{
        exgcd(b,a%b,d,y,x);
        y-=x*(a/b);
    }
}
ull inv(ull a){
    ull d,x,y;
    exgcd(a,R,d,x,y);
    return (x+R)%R;
}
int main(){
    int T;
    int n,m;
    scanf("%d%llu",&T,&R);
    sievePrime();
    sieveFac();
    sievePhiFac();
    while(T--){
        scanf("%d%d",&n,&m);
        printf("%llu\n",(phifac[m]*((fac[n]*inv(fac[m]))%R))%R);
    }
    return 0;
}