每次给出 \(n,k,p\),求出长为 \(n\) 的回文串以及其旋转变换的总数,且字符集大小为 \(k\),答案对 \(p\) 取模。
\(T\le 10\),\(n\le 10^{18}\),\(k\le n\),\(10^9\le p\le 2^{30}\).
首先长为 \(n\) 字符集大小为 \(k\) 的回文串总数显然是 \(\displaystyle g(n)=k^{\lfloor\frac{n+1}{2}\rfloor}\).
然后思考这些串有多少不同的旋转变换,这个时候就要枚举最小循环节。
对于循环节长度 \(d\),若 \(d\) 为奇数,串的贡献就是 \(d\).
然后看下偶数的情况,容易发现 abcd abcd \(\rightarrow\) cdab cdab 的构造,也就是说一个串的贡献应该是 \(\displaystyle\frac{d}{2}\) 了。
把这个函数刻画成 \(\displaystyle h(d)=\begin{cases}d&d\operatorname{mod}2=1\\\frac{d}{2}&d\operatorname{mod}2=0\end{cases}\)
设最小循环节长度为 \(d\) 的回文串数量为 \(f(d)\),枚举 \(d\) 得答案即
\[\sum_{d|n}h(d)f(d) \]问题是还要求这个 \(f\).
有一个比较显然的性质:
\[g(m)=k^{\lfloor\frac{m+1}{2}\rfloor}=\sum_{d|m}f(d) \]\[g=f * 1\rightarrow f=g * \mu \]原式即
\[\sum_{d|n}h(d)\sum_{p|d}g(p)\mu(\frac{d}{p}) \]\[\sum_{p|n}g(p)\sum_{d|\frac{n}{p}}h(dp)\mu(d) \]发现必须得把 \(h(dp)\) 给拆开,想一下何时有 \(h(dp)=dh(p)\).
上式不成立当且仅当 \(d\) 为偶数且 \(p\) 为奇数,那么 \(\displaystyle\frac{n}{p}\) 也是偶数。
故 \(\displaystyle\sum_{d|\frac{n}{p}}h(dp)\mu(d)=p\sum_{d|\frac{n}{p}}h(d)\mu(d)\).
考虑有贡献的位,即 \(\mu(d)\) 非 \(0\) 的 \(d\).
因为 \(\displaystyle\frac{n}{p}\) 是偶数,那么奇数 \(d\) 一定对应其 \(2d\),且有 \(h(k)=h(2k)\),\(\mu(d)=-\mu(2d)\),这两项加起来就抵消了,所以此时求和式的值为 \(0\).
先把这一 part 的贡献忽略掉再看看。
\[\sum_{p|n}g(p)h(p)\sum_{d|\frac{n}{p}}d\mu(d) \]再看一下怎么比较快地求 \(\displaystyle\sum_{d|n}d\mu(d)\).
容易发现这个是有组合意义的,只需要考虑所有互不相同的质因子。将 \(p\) 加入时,选 \(p\) 则当前贡献乘上 \(-p\),否则不变,那么 \(f=(f'-pf')=(1-p)f'\),可以得到
\[\sum_{d|n}d\mu(d)=\prod_{p\in\mathbb P,p|n}(1-p) \]这个式子是枚举 \(p\) 时补集上的东西。
总结一下式子
\[\sum_{d|n}g(d)h(d)\prod_{p\in\mathbb P,p|\frac{n}{d}}(1-p) \]由于 \(\displaystyle\max_{d=1}^{\rm 1e18}\{\sigma_0(d)\}=103680\),用 Pollard Rho 分解 \(n\),搞出来分式上面的东西,就可以枚举质因子跑 dfs 了。
后面一半式子要提出来做。
记得要跳过 \(x\) 为奇数且 \(\displaystyle\frac{n}{x}\) 为偶数的情况。
默认三者同阶,单次时间复杂度 \(O(n^{\frac{1}{4}}\log n+\sigma_0(n)\log n)\).
时间上面相对来说还是很友好的。这种玩意还是第一次手写。
这种能够对奇偶讨论的东西实际上真的很好用。
#include<bits/stdc++.h>
#define ll long long
#define Test 10
#define Omega 16
#define Sigma 103681
using namespace std;
mt19937 rd(233);
ll read(){
ll x=0,w=1;char ch=getchar();
while(!isdigit(ch)){if(ch=='-')w=-1;ch=getchar();}
while(isdigit(ch))x=(x<<3)+(x<<1)+(ch^48),ch=getchar();
return x*w;
}
ll qmul(ll a,ll b,ll p){
ll r=((long double)a/p*b+0.5);
r=a*b-r*p;
return r<0?r+p:r;
}
ll qpow(ll k,ll b,ll p){
ll ret=1;k%=p;
while(b){
if(b&1)ret=qmul(ret,k,p);
k=qmul(k,k,p),b>>=1;
}
return ret;
}
bool check(ll a,ll n){
ll u=n-1,t=0;
while(!(u&1))u>>=1,t++;
ll x=qpow(a,u,n);
if(x==1)return false;
for(ll i=1,tp;i<=t;i++,x=tp){
tp=qmul(x,x,n);
if(x!=1&&x!=n-1&&tp==1)
return true;
}
return x!=1;
}
bool Miller_Rabin(ll n){
if(n==2)return true;
if(n<2||!(n&1))return false;
for(int i=1;i<=Test;i++){
if(check(rd()%(n-1)+1,n))
return false;
}
return true;
}
ll f(ll x,ll c,ll n){
return (qmul(x,x,n)+c)%n;
}
ll Pollard_Rho(ll n){
if(!(n&1))return 2;
ll x=rd()%(n-1)+1,y=x,c=rd()%(n-1)+1,ret=1;
for(ll i=1;;i<<=1,y=x,ret=1){
for(ll j=1;j<=i;j++){
x=f(x,c,n);
ret=qmul(ret,abs(x-y),n);
if(!(j&127)){
ll d=__gcd(ret,n);
if(d>1)return d;
}
}
ll d=__gcd(ret,n);
if(d>1)return d;
}
}
ll pr[Omega];int cnt[Omega],tot;
map<ll,int>mp;
void find(ll n,int now){
if(n==1)return;
if(Miller_Rabin(n)){
if(!mp[n]){
mp[n]=++tot,pr[tot]=n;
}
cnt[mp[n]]+=now;
return;
}
ll p=n;
while(p>=n)p=Pollard_Rho(n);
int cur=0;
while(n%p==0)n/=p,cur++;
find(p,now*cur),find(n,now);
}
ll n,k,P,ans;
ll g(ll d){
return qpow(k,(d+1)>>1,P);
}
ll h(ll d){
return (d&1)?d%P:(d>>1)%P;
}
ll _M[Sigma],Sum;
map<ll,ll>dv;
void dfs(int id,int cur,ll M,ll S){
if(id==tot+1){
_M[++Sum]=M,dv[M]=S;
return;
}
if(cur<cnt[id])dfs(id,cur+1,M*pr[id],S*(!cur?(1-pr[id]):1));
dfs(id+1,0,M,S);
}
int main(){
int T=read();
while(T--){
for(int i=1;i<=tot;i++)pr[i]=cnt[i]=0;
mp.clear(),tot=ans=Sum=0,dv.clear();
n=read(),k=read(),P=read();
find(n,1);
dfs(1,0,1,1);ll M;
for(int i=1;i<=Sum;i++){
M=_M[i];
if((M&1)&&(!((n/M)&1)))continue;
(ans+=(dv[n/M]%P+P)%P*g(M)%P*h(M)%P)%=P;
}
printf("%lld\n",ans);
}
return 0;
}
不知道为什么山东一年两道困难反演题。
标签:frac,P4607,ll,displaystyle,mu,SDOI2018,return,sum,回文 From: https://www.cnblogs.com/SError0819/p/17621537.html