[bzoj3202][SDOI2013]项链

3202: [Sdoi2013]项链

Time Limit: 30 Sec Memory Limit: 512 MB
Submit: 392 Solved: 135
[Submit][Status][Discuss]
Description

项链是人体的装饰品之一，是最早出现的首饰。项链除了具有装饰功能之外，有些项 链还具有特殊显示作用，如天主教徒的十字架
链和佛教徒的念珠。 从古至今人们为了美化人体本身，也美 化环境，制造了各种不同风格，不同特点、不同式样的项链，满足了不同肤色、不同民族、不同审美观的人的审美需要。就材料而论，首
饰市场上的项链有黄金、白银、珠宝等几种。珍珠项链为珍珠制成的饰品，即将珍珠 钻孔后用线串在一起，佩戴于项间。天然珍珠项链具有一定的护养作用。

最近，铭铭迷恋上了一种项链。与其他珍珠项链基本上相同，不过这种项链的珠子却 与众不同，是正三菱柱的泰山石雕刻而成的。三菱柱的侧面是正方形构成的，上面刻有数字。 能够让铭铭满意的项链必须满足下面的条件：
1：这串项链由n颗珠子构成的。
2：每一个珠子上面的数字x,必须满足0 < x <= a，且珠子上面的数字的最大公约数要恰 好为1。两个珠子被认为是相同的，当且仅当他们经过旋转，或者翻转后能够变成一样的。 3：相邻的两个珠子必须不同。
4：两串项链如果能够经过旋转变成一样的，那么这两串项链就是相同的！ 铭铭很好奇如果给定n和a，能够找到多少不同串项链。由于答案可能很大，所以对输 出的答案mod 1000000007。

Input

数据由多组数据构成： 
第一行给定一个T<=10，代表由T组数据。 
接下来T行，每行两个数n和a。 

Output

对于每组数据输出有多少不同的串。 

Sample Input

1 

2  2

Sample Output

3 

HINT

对于100%的数据：所有的n<=10^14，a<=10^7，T<=10； 

样例解释：由三种珠子：[1,1,1],[1,1,2],[1,2,2].组成的串有：[1,2],[1,3],[2,3]。

先不考虑题目中所说的$gcd=1$的情况，那么就是从可重元素中取k个有多少种取法。因为可以重复选，所以多加入$k$个元素，就相当于要把这$n+k$个元素分成$n$份有多少种取法。可以用插板法，在$n+k-1$个空中找$n-1$个空。就是$C(n+k-1,n-1)$
也就等价于我们第一问要求的珠子的种数$o=C(n+2,3)$
然后再把$gcd!=1$的充斥掉：那就是减去1个质数的，加上2个质数的……，那么最后就是用莫比乌斯函数容斥一下：$o=\sum_{d|a}\mu(d)*C(\frac{a}{d}+2,3)$
然后就应该用到群论的知识了，由于不会群论，所以听了大神莫比乌斯反演的方法。
$F(x)表示最短循环节长度为x的数量$
$f(x)=\sum_{i|x}F(x)$
那么最后的答案就是$ans=\sum_{d|n}\frac{F(d)}{d}$因为旋转之后相同的算一种，所以一个长度为d的会被统计d次。
$ans=\sum_{d|n}\frac{F(d)}{d}$
$\ \ \ \ \ \ =\sum_{d|n}\frac{\sum_{p|d}{\mu(p)f(\frac{d}{p})}}{d}$
$\ \ \ \ \ \ =\sum_{d|n}\frac{1}{d}\sum_{p|d}\mu(p)f(\frac{d}{p})$
$\ \ \ \ \ \ =\sum_{t|n}\sum_{p|\frac{n}{t}}\frac{1}{pt}\mu(p)f(t)$
$\ \ \ \ \ \ =\sum_{t|n}\frac{f(t)}{t}\sum_{p|\frac{n}{t}}\frac{\mu(p)}{p}$
$根据容斥原理：\varphi(n)=\sum_{d|n}\frac{\mu(d)*n}{d}$
$所以\sum_{p|\frac{n}{t}}\frac{\mu(p)}{p}=\frac{\varphi(\frac{n}{t})}{\frac{n}{t}}$
$ans=\sum_{t|n}\frac{f(t)}{t}\frac{\varphi(\frac{n}{t})}{\frac{n}{t}}$
$\ \ \ \ \ \ =\frac{1}{n}\sum_{t|n}f(t)\varphi(\frac{n}{t})$
现在我们只需要求出f(t)就行了。
设g[i][0]表示长度为i的开头和末尾相同的序列的个数，g[i][1]表示表示长度为i的开头和末尾不相同的序列的个数。
那么最后的答案$f(n)=g[n][1]$。f表示的是一个和的形式，但是在求g[n][1]的过程中就已经统计过了长度是n的约数的序列的个数了。
比如在统计长度为6的时候，就会把$1 2 3 1 2 3$这种序列统计进去了。素以这样求是正确的。
$g[i][0]=g[i-1][1]$
$g[i][1]=(o-2)g[i][1]+(o-1)g[i-1][0]$
$\ \ \ \ \ \ \ \ \ =(o-2)g[i][1]+(o-1)g[i-2][1]$
$所以f(n)=(o-2)f(n-1)+(o-1)f(n-2)$
这个东西可以矩乘，也可以用特征根求通项。
$r^2=(o-2)r+(o-1)$
$解出来r1=-1,r2=d-1$
然后$f(n)=c1*r1^n+c2*r2^n$
可以把$f(0)=0,f(1)=n*(n-1)$带入就可以求出$c1,c2$
$解得c1=d-1,c2=1$
$f(n)=(d-1)(-1)^n+(d-1)^n$
直接dfs求出所有的约数算就可以了。
但是这个题vfk加强了数据！！！可能会有n是p的倍数的情况，也就是n没有逆元。怎么办呢?
因为$\frac{x}{y}\%p=\frac{x\%yp}{y}$
我们看要求的答案是$\frac{\sum_{t|n}f(t)\varphi(\frac{n}{t})}{n}\%p$用x代替$\sum$那么就是
$ans=\frac{x}{n}\%p=\frac{\frac{x}{p}\%p^2}{\frac{n}{p}}\%p$

#include<iostream>
#include<cstdio>
#include<cstring>
using namespace std;
#define O 1000000007
#define LL long long
const int N=10000000;
bool flag[N+10];
LL o,n,ni,phi,six,D=1e9+7,p[N],ans,use[N];
int prime[N+10],u[N+10],T,m,mi[N];
inline void prepare(){
    int i,j;
    for(u[1]=1,i=2;i<=N;++i){
        if(!flag[i]){
            prime[++prime[0]]=i;
            u[i]=-1;
        }
        for(j=1;j<=prime[0]&&i*prime[j]<=N;++j){
            flag[i*prime[j]]=true;
            if(i%prime[j]==0){
                u[i*prime[j]]=0;
                break;
            }
            u[i*prime[j]]=-u[i];
        }
    }
    for(i=2;i<=N;++i) u[i]+=u[i-1];
}
inline LL quickmul(LL x,LL y){
    LL sum=0,v=1;
    if(y<0LL) y=-y,v=-1;
    while(y){
        if(y&1LL) sum=(sum+x)%D;
        y>>=1LL;x=(x+x)%D;
    }
    return sum*v;
}
inline LL quickpow(LL x,LL y){
    LL sum=1LL;
    while(y){
        if(y&1L) sum=quickmul(sum,x);
        y>>=1LL;x=quickmul(x,x);
    }
    return sum;
}
inline LL C(int x){
    return quickmul(quickmul(x,x-1),quickmul(x-2,six));
}
inline LL calc(){
    int i,j;
    LL sum=1,f,PHI=1;
    for(i=1;i<=p[0];++i)
      for(j=1;j<=use[i];++j) sum=sum*p[i];
    for(PHI=sum,i=1;i<=p[0];++i)
      if(use[i]) PHI=PHI/p[i]*(p[i]-1);
    sum=n/sum;
    f=(((o-1LL)*(LL)(sum%2?-1:1)+quickpow((LL)(o-1),sum))+D)%D;
    return quickmul(PHI,f);
}
inline void dfs(int x){
    int i;
    if(x==p[0]+1){
        ans=(ans+calc())%D;
        return ;
    }
    for(i=0;i<=mi[x];++i)
      use[x]=i,dfs(x+1);
}
int main(){
    prepare();
    scanf("%d",&T);
    while(T--){
        LL pre=O;D=O;
        int i,pos,mark=0;
        scanf("%lld%d",&n,&m);
        if(n%D) phi=D-1,ni=quickpow(n,D-2);
        else ni=quickpow(n/D,D-2),D=D*D,phi=quickmul(pre,pre-1),mark=1;
        six=quickpow(6LL,phi-1);
        for(o=0,i=1;i<=m;i=pos+1){
            pos=m/(m/i);
            o=(o+(quickmul((LL)(u[pos]-u[i-1]),C(m/i+2))+D)%D)%D;
        }
        LL x=n;
        for(p[0]=mi[0]=0,i=1;i<=prime[0]&&x>1LL;++i){
            if(x%prime[i]==0){
                p[++p[0]]=prime[i];
                mi[++mi[0]]=0;
                while(x%prime[i]==0){
                    ++mi[mi[0]];
                    x/=prime[i];
                }
            }
        }
        if(x>1LL) p[++p[0]]=x,mi[++mi[0]]=1;
        ans=0;dfs(1);
        if(mark) ans/=pre;
        ans=quickmul(ans,ni);
        printf("%lld\n",ans%pre);
    }
}