luogu
bzoj
sol
枚举一个长度为\(n\)为回文串，它的所有循环位移都可以产生贡献。
但是这样算重了。重复的地方在于可能多个回文串循环同构，或者可能有的回文串经过小于\(n\)次循环位移后能够得到自身。
一个比较好的处理方式是：对每个回文串求最小的\(x\)使这个串经过\(x\)次循环位移后可以再次成为一个回文串。这样对每个回文串求\(\sum x\)显然就不会算重了。
考虑一个串的\(x\)是什么。显然会和这个串的最小循环节长度有关。实际上如果最小循环节长度为偶数，那么\(x\)就会是这个长度的一半；否则就等于这个长度。
形式化地，如果一个回文串的最小循环节长度为\(i\)，那么它对答案的贡献就是\(h(i)=i\frac{1+[i\mbox{是奇数}]}{2}\)。
设最小循环节为\(i\)的回文串共有\(f(i)\)个，那么我们要求的答案就是
\[Ans=\sum_{d|n}f(d)h(d)\]
又因为\[\sum_{d|n}f(d)=k^{\lceil\frac n2\rceil}=g(n)\]
所以\[f(n)=\sum_{d|n}g(d)\mu(\frac nd)\]
代入原式\[Ans=\sum_{d|n}\sum_{i|d}g(i)\mu(\frac di)h(d)\\=\sum_{i|n}g(i)\sum_{d|\frac ni}\mu(d)h(id)\]
我们希望可以把\(h(id)\)中的\(i\)提出来，这样后半部分就是一个关于\(\frac ni\)的函数了。
因为\(h(x)\)不是\(x\)就是\(\frac x2\)，我们发现\(h(id)\neq d\times h(i)\)当且仅当\(i\)是奇数且\(d\)是偶数，而\(d|\frac ni\)所以\(d\)是偶数就说明\(\frac ni\)也是偶数。那么我们现在假设\(i\)是奇数且\(\frac ni\)是偶数，考虑下面这个式子的取值。
\[\sum_{d|\frac ni}\mu(d)h(id)\]
显然只有\(\mu(d)\)非零项有贡献，而\(\frac ni\)中含有\(2\)这个因子就使得所有\(\mu(d)\)非零项中含\(2\)与不含\(2\)的\(d\)可以一一对应。他们的\(h(id)\)的值是相同的，而\(\mu(d)\)的值恰好相反，所以这个式子的值一定为\(0\)。
话说回来。我们现在已经知道了\(h(id)\neq d\times h(i)\)的情况没有贡献，就可以放心大胆地用这一种变换了。
\[Ans=\sum_{i|n}g(i)\sum_{d|\frac ni}\mu(d)h(id)\\=\sum_{i|n}g(i)h(i)\sum_{d|\frac ni}d\mu(d)\]
（在枚举\(i\)时需要跳过\(i\)是奇数而\(\frac ni\)是偶数的项）
考虑后面的东西是个啥。还是只有\(\mu(d)\)非零项有贡献，也就是说含有奇数个质因子的\(d\)会乘上\(-1\)，含有偶数个质因子的\(d\)为乘上\(1\)，所以这个值相当于是将\(\frac ni\)质因数分解为\(p_1^{a_1}p_2^{a_2}...p_k^{a_k}\)后，为\(\prod_{i=1}^k(1-p_i)\)。
所以到这里就比较简单了。先用\(Pollard-Rho\)算法将\(n\)分解，再dfs枚举\(n\)的每一个约数\(d\)，在搜索的过程中自然可以求出那个\(\prod_{i=1}^k(1-p_i)\)。
code
#include<cstdio>
#include<algorithm>
#include<ctime>
using namespace std;
#define ll long long
ll mul(ll x,ll y,ll m){x%=m;y%=m;return (x*y-(ll)(((long double)x*y+0.5)/(long double)m)*m+m)%m;
}
ll fastpow(ll x,ll y,ll m){ll res=1;while (y) {if (y&1) res=mul(res,x,m);x=mul(x,x,m);y>>=1;}return res;
}
ll f[]={2,3,5,7,11,13,17,19,23,29};
bool MR(ll p){for (int i=0;i<10;++i){if (p<=f[i]) break;if (fastpow(f[i],p-1,p)!=1) return false;ll pp=p-1;while (~pp&1){pp>>=1;ll y=fastpow(f[i],pp,p);if (mul(y,y,p)==1&&y!=1&&y!=p-1) return false;}}return true;
}
ll PR(ll n,ll c){ll i=0,k=2,x,y;x=y=1+rand()%(n-1);while (1){x=(mul(x,x,n)+c)%n;ll d=__gcd((y-x+n)%n,n);if (d!=1&&d!=n) return d;if (x==y) return n;if (++i==k) y=x,k<<=1;}
}
ll tmp[100];int len;
void fact(ll n){if (n==1) return;if (MR(n)) {tmp[++len]=n;return;}ll p=n;for (int c=233;p==n;--c) p=PR(p,c);fact(p);fact(n/p);
}
int Case,q[100],cnt,mod;ll n,k,p[100],ans;
int fpow(int x,ll y){ll res=1;while (y) {if (y&1) res=1ll*res*x%mod;x=1ll*x*x%mod;y>>=1;}return res;
}
int g(ll n){return fpow(k,(n+1)>>1);}
int h(ll n){return (n&1?n:n>>1)%mod;}
void dfs(int i,ll d,int pro){if (i==cnt+1){if ((n/d&1)&&(d&1)==0) return;ans=(ans+1ll*g(n/d)*h(n/d)%mod*pro)%mod;return;}dfs(i+1,d,pro);pro=1ll*pro*(mod+1-p[i]%mod)%mod;for (int j=1;j<=q[i];++j) d*=p[i],dfs(i+1,d,pro);
}
int main(){srand(20020415);scanf("%d",&Case);while (Case--){scanf("%lld%lld%d",&n,&k,&mod);k%=mod;len=cnt=0;fact(n);sort(tmp+1,tmp+len+1);for (int i=1;i<=len;++i){if (tmp[i]!=tmp[i-1]) p[++cnt]=tmp[i],q[cnt]=0;++q[cnt];}ans=0;dfs(1,1,1);printf("%lld\n",ans);}return 0;
}
