题目链接:https://www.luogu.com.cn/problem/P3307
n个珠子的一个环形项链,每个珠子有三个1∼k的整数。
- 两个珠子不同当且仅当它们不能通过翻转或者旋转得到
- 两个项链不同当且仅当它们不能通过旋转得到
- 珠子要求上面的数字互质
- 项链要求相邻珠子不同
求方案数。
珠子的计数和项链计数是没什么关系的,所以两个分开求。
先考虑求珠子有多少种。相当与求有多少种三元组,因为三个数字的排列通过以上操作可以得到任何其他排列。
首先答案包括了三个不同的,有两个相同的,全部都相同的。
也就是设Si表示i元组的数目的话答案就是S3+3S2+2S1(容斥一下可以得到的)
可以先不考虑顺序问题,最后除上一个六就好了。
要求互质所以可以直接上莫反就变为求∑ki=1μ(i)F(i)。对于每个F(i)来说就是有⌊ki⌋个数字不要求互质的答案。
对于三元组来说F(i)=⌊ki⌋3,二元组就是F(i)=⌊ki⌋2可以和三元组一起解决,一元组就是1。
这样整除分块就可以√k的处理单次询问了。
现在考虑项链的方案,根据Burnside引理的话就是1n∑ni=1f(gcd(i,n)),其中f(i)表示大小为i的不动点数量(也就是转i次的话置换中环的数量)。化一下就是1n∑d|nf(d)φ(nd),需要考虑如何求f。
因为是不动点,所以每一个环中的颜色都得相同,然后因为每个环都是相互隔开的,也就是可以视为有gcd(n,i)个珠子头尾相连要涂色,不能有相邻颜色相同。
这个f是有递推式f(i)=f(i−1)∗(m−2)+f(i−2)∗(m−1)(据说法是第一个是直接在末端插入一个颜色的珠子,第二个是先在末端插入一个和首部相同颜色的珠子再插入一个不同颜色的将它们隔开,这个是第一个无法处理的情况)。
这个可以用矩阵乘法快速递推,但是好像有更简单的方法可以一试,考虑将式子化为f(x)−af(x−1)=bf(x−1)+abf(x−2)的形式,这样可以将答案化为一个带幂的函数。
那就有方程{b−a=m−2ab=m−1,然后可以解出a=1或者a=1−m。带个1进去做。也就是有f(x)−f(x−1)=(m−1)(f(x−1)+f(x−2))也就是乘上m−1就往后退一步,所以有f(x)−f(x−1)=(m−1)x−2(f(1)+f(2))直接带f(1)=0,f(2)=m(m−1)进去就有f(x)=m(m−1)x−1−f(x−1)
设一个A满足
f(x)+(m−1)A=−(f(x−1)+A)
f(x)+mA=−f(x−1)
又因为有f(x)=m(m−1)x−1−f(x−1)
m(m−1)x−1+mA=0⇒A=(m−1)x−1
所以把A带回最初式子,把f(x−1)换成f(x)就有式子
f(x)+(m−1)x=−(f(x)−(m−1)x−1)
⇒f(x)=(−1)x(m−1)+(m−1)x
(这个好像是叫特征方程的玩意)
现在就可以用快速幂算f(x)了,φ比较大,不能直接√n求单个,可以先把n质因数分解了然后dfs下去就好了。然后还有因为n很大,有可能是模数P的倍数,这个时候我们求答案模P2,然后除以一个P在和nP搞逆元就好了。
时间复杂度O(k+T(√n+√m))
#include<cstdio>
#include<cstring>
#include<algorithm>
#define ll long long
using namespace std;
const ll N=1e7+10,XJQ=1e9+7,inv6=833333345000000041ll;
ll T,n,m,k,ans,P,cnt,num;
ll pri[N/10],mu[N],p[110],c[110];
bool v[N];
void init(){
mu[1]=1;
for(ll i=2;i<N;i++){
if(!v[i])pri[++num]=i,mu[i]=-1;
for(ll j=1;j<=num&&i*pri[j]<N;j++){
v[i*pri[j]]=1;
dans
}
for(ll i=1;i<N;i++)
mu[i]+=mu[i-1];
return;
}
ll ksc(ll a,ll b,ll p){
a%=p;b%=p;
ll tmp=(long double)a*b/p;
long double ans=a*b-tmp*p;
if(ans>=p)ans-=p;
else if(ans<0)ans+=p;
return ans;
}
ll power(ll x,ll b,ll p){
ll ans=1;
while(b){
if(b&1)ans=ksc(ans,x,p);
x=ksc(x,x,p);b>>=1;
}
return ans;
}
void Get_M(ll n){
m=2;
for(ll l=1,r;l<=n;l=r+1){
r=n/(n/l);ll k=n/l;
(m+=ksc(ksc(ksc(k,k,P),k+3,P),mu[r]-mu[l-1],P))%=P;
}
m=ksc(m,inv6,P);
return;
}
ll f(ll n)
{return power(m-1,n,P)+((n&1)?(P-m+1):(m-1))%P;}
void dfs(ll dep,ll x,ll val){
if(dep>cnt){
(ans+=ksc(f(n/x),val,P))%=P;
return;
}
dfs(dep+1,x,val);
val*=p[dep]-1;x*=p[dep];
dfs(dep+1,x,val);
for(ll i=2;i<=c[dep];i++)
val=val*p[dep],x*=p[dep],dfs(dep+1,x,val);
return;
}
signed main()
{
scanf("%lld",&T);
init();
while(T--){
scanf("%lld%lld",&n,&k);
if(n%XJQ!=0)P=XJQ;
else P=XJQ*XJQ;
Get_M(k);
ll x=n;cnt=ans=0;
for(ll i=1;i<=num&&pri[i]*pri[i]<=x;i++){
if(x%pri[i]==0){
p[++cnt]=pri[i];c[cnt]=0;
while(x%pri[i]==0)
c[cnt]++,x/=pri[i];
}
}
if(x>1)p[++cnt]=x,c[cnt]=1;
dfs(1,1,1);
if(n%XJQ==0)
printf("%lld\n",ans/XJQ*power(n/XJQ,XJQ-2,XJQ)%XJQ);
else
printf("%lld\n",ans*power(n,XJQ-2,XJQ)%XJQ);
}
return 0;
}
__EOF__
【推荐】国内首个AI IDE,深度理解中文开发场景,立即下载体验Trae
【推荐】编程新体验,更懂你的AI,立即体验豆包MarsCode编程助手
【推荐】抖音旗下AI助手豆包,你的智能百科全书,全免费不限次数
【推荐】轻量又高性能的 SSH 工具 IShell:AI 加持,快人一步
· .NET Core 中如何实现缓存的预热?
· 从 HTTP 原因短语缺失研究 HTTP/2 和 HTTP/3 的设计差异
· AI与.NET技术实操系列:向量存储与相似性搜索在 .NET 中的实现
· 基于Microsoft.Extensions.AI核心库实现RAG应用
· Linux系列:如何用heaptrack跟踪.NET程序的非托管内存泄露
· TypeScript + Deepseek 打造卜卦网站:技术与玄学的结合
· 阿里巴巴 QwQ-32B真的超越了 DeepSeek R-1吗?
· 【译】Visual Studio 中新的强大生产力特性
· 张高兴的大模型开发实战:(一)使用 Selenium 进行网页爬虫
· 【设计模式】告别冗长if-else语句:使用策略模式优化代码结构