[BZOJ3202][SDOI2013]项链

bzoj
luogu

sol

数论多合一。
可以考虑分两步做。
1、求有多少种不同的珠子,设其数量为\(m\)
2、用\(m\)种颜色对\(n\)个珠子的项链染色,要求相邻不同色,且旋转相同视作同一种方案,求方案数。

求有多少种不同的珠子

\(A=\sum_{i=1}^a\sum_{j=1}^a\sum_{k=1}^a[\gcd(i,j,k)=1]\)
然而算重了。
如果\(i,j,k\)都不相等的话就多算了\(6\)次,如果有两个相等就多算了\(3\)次,如果三个都相等......那肯定是三个\(1\)啦,并没有多算。
这个问题怎么解决呢?我们再算\(B=\sum_{i=1}^a\sum_{j=1}^a[\gcd(i,j)=1]\),那么答案就是\(\frac 16(A+3B+2)\)(把每种情况补齐到\(6\)次,再一起除以\(6\))。
\(A\)\(B\)是很简单的莫比乌斯反演,具体来说就是:
\[A=\sum_{i=1}^a\lfloor\frac ai\rfloor^3\mu(i)\\B=\sum_{i=1}^a\lfloor\frac ai\rfloor^2\mu(i)\]

求染色方案

这里要用到\(Burnside\)定理。
这里的置换只有一种,就是旋转。设旋转\(i\)个位置,那么就会有产生\(\gcd(i,n)\)个等价类。
所以答案就是\[\frac{\sum_{i=1}^nf_{\gcd(i,n)}}{n}\]其中\(f_i\)表示\(i\)个点的环染色的方案数。也可以认为是给一个长为\(n\)的序列染色,然后保证首尾不同色的方案数。
考虑对上式做一些简化。因为\(\gcd(i,n)\)一定是\(n\)的约数,所以可以枚举某个约数\(d\),然后答案就变成了\[\frac{\sum_{d|n}f_d\sum_{i=1}^n[\gcd(i,n)=d]}{n}=\frac{\sum_{d|n}f_d\varphi(\frac nd)}{n}\]

\(f_n\)

考虑一下\(f_n\)的递推。
\[f_n=(m-2)f_{n-1}+(m-1)f_{n-2},f_1=0,f_2=m(m-1)\]
讲道理的话,前一种转移表示在末尾填一个与首尾都不相同的颜色,这样有\(m-2\)种选择;后一种转移表示先填一个和首位相同的颜色,再填一个不同色,这样有\(1\times(m-1)\)种选择。
然后就可以矩乘了。

\(f_n\)的通项公式

有人跟我讲这题要求通项公式我就没去写矩乘。常数优秀一点说不定跑得过?
下面讲一下求通项公式
\[f_n-(m-2)f_{n-1}-(m-1)f_{n-2}=0\]
其特征方程\[x^2-(m-2)x-(m-1)=0\]
解得其特征根\[x_1=-1,x_2=m-1\]
于是设其通项公式为\[f_i=a(-1)^i+b(m-1)^i\]
\(f_1=0,f_2=m(m-1)\)代入解得\[a=1-m,b=1\]
所以\[f_i=(1-m)(-1)^i+(m-1)^i\]
复杂度还是\(O(\log n)\)但是常数小很多。

\(\varphi(n)\)

最好不要枚举因数然后暴力算\(\varphi(n)\)
一种可行的方案是:先将\(n\)质因数分解,然后\(dfs\)搜出每一个因数,在搜的过程中自然可以求出它的\(\varphi\)
这样的话总的复杂度就是\(O(a+T(\sqrt a+d(n)\log n))\)

你以为这样就做完了吗

\(n\le10^{14}\)
这样可能会产生一个问题:\(n\)\(P=10^9+7\)的倍数。这样算出来的分子要除以\(n\),不好意思不存在逆元。
怎么办呢?
我们把答案对\(P^2\)取模。如果\(P|n\),那么分子既然是\(n\)的倍数就也是\(P\)的倍数。直接除以\(P\),接着还需要除以\(\frac nP\),而\(\frac nP\)在模\(P\)意义下是存在逆元的(\(n\le P^2\)),于是就可以乘逆元了。

code

#include<cstdio>
#include<algorithm>
using namespace std;
#define ll long long
const ll mod = 1000000007ll;
const ll MOD = 1000000007ll*1000000007ll;
const ll Inv6 = 833333345000000041ll;
const int N = 10000005;
int T,zhi[N],pri[N/10],tot,mu[N],cnt,q[30];ll n,m,p[30],ans;
void add(ll &x,ll y){x+=y;if(x>=MOD)x-=MOD;}
ll mul(ll x,ll y){
    return (x*y-(ll)(((long double)x*y+0.5)/(long double)MOD)*MOD+MOD)%MOD;
}
ll sqr(ll x){return mul(x,x);}
ll cub(ll x){return mul(mul(x,x),x);}
ll fastpow(ll a,ll b){
    ll res=1;
    while (b) {if (b&1) res=mul(res,a);a=mul(a,a);b>>=1;}
    return res;
}
ll fpow(ll a,ll b){
    ll res=1;
    while (b) {if (b&1) res=res*a%mod;a=a*a%mod;b>>=1;}
    return res;
}
void init(){
    mu[1]=1;
    for (int i=2;i<N;++i){
        if (!zhi[i]) pri[++tot]=i,mu[i]=-1;
        for (int j=1;i*pri[j]<N;++j){
            zhi[i*pri[j]]=1;
            if (i%pri[j]) mu[i*pri[j]]=-mu[i];
            else break;
        }
    }
    for (int i=1;i<N;++i) mu[i]+=mu[i-1];
}
ll cal(ll n){
    ll A=0,B=0;
    for (ll i=1,j;i<=n;i=j+1){
        j=n/(n/i);
        add(A,mul(cub(n/i),mu[j]-mu[i-1]+MOD)%MOD);
        add(B,mul(sqr(n/i),mu[j]-mu[i-1]+MOD)%MOD);
    }
    add(A,mul(B,3));add(A,2);return mul(A,Inv6);
}
ll F(ll n){
    ll res=n&1?1-m:m-1;add(res,MOD);
    add(res,fastpow(m-1+MOD,n));return res;
}
void dfs(int i,ll cur,ll phi){
    if (i==cnt+1) {add(ans,mul(phi,F(n/cur)));return;}
    dfs(i+1,cur,phi);
    cur*=p[i];phi*=p[i]-1;dfs(i+1,cur,phi);
    for (int j=2;j<=q[i];++j)
        cur*=p[i],phi*=p[i],dfs(i+1,cur,phi);
}
int main(){
    init();
    scanf("%d",&T);while (T--){
        scanf("%lld%lld",&n,&m);m=cal(m);
        ll x=n;cnt=0;
        for (int i=1;i<=tot&&(ll)pri[i]*pri[i]<=x;++i)
            if (x%pri[i]==0){
                p[++cnt]=pri[i];q[cnt]=0;
                while (x%pri[i]==0) x/=pri[i],++q[cnt];
            }
        if (x>1) p[++cnt]=x,q[cnt]=1;
        ans=0;dfs(1,1,1);
        if (n%mod==0){
            ans/=mod;
            ans=ans*fpow(n/mod,mod-2)%mod;
        }else{
            ans%=mod;
            ans=ans*fpow(n%mod,mod-2)%mod;
        }
        printf("%lld\n",ans);
    }
    return 0;
}

转载于:https://www.cnblogs.com/zhoushuyu/p/9657640.html

  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值