HDU 6427 Problem B. Beads（polya+数论+素因子分解Pollard-rho）

Description

用$m$种颜色给一串长度为$n$的项链染色，旋转和翻转视作一种方案，且颜色平移也视为一种方案，问染色方案数

Input

第一行一整数$T$表示用例组数，每组用例输入两个整数$n,m$

$(1\le T\le 30,3\le n\le 10^{18},2\le m\le 10^{18},n,m\not|998244353)$

Output

输出染色方案数，结果模$998244353$

Sample Input

5
3 2
4 2
8 5
9 5
2333 333

Sample Output

2
4
5079
22017
544780894

Solution

$polya$，变换群大小为$2nm$，对于旋转和颜色平移的$nm$个变换，若一个染色方案能够在颜色平移为$g$（即$1$颜色变成$g$颜色），旋转平移为$d$（即$1$位置移动到$d$位置）后不变，那么从$1$走回$1$的这$\frac{n}{gcd(d,n)}$个位置，其颜色必然要整数次从$1$变成$1$，而$1$颜色变回$1$颜色至少需要$\frac{m}{gcd(g,m)}$次，故这$nm$个变换下不变的染色方案数为

$$\begin{array}{lcl} \sum\limits_{d=1}^nm^{gcd(d,n)}\sum\limits_{g=1}^m[\frac{m}{gcd(g,m)}|\frac{n}{gcd(d,n)}]&=&\sum\limits_{d|n}\varphi(d)m^{\frac{n}{d}}\sum\limits_{g|m}\varphi(g)[g|d]\\ &=&\sum\limits_{d|n}\varphi(d)m^{\frac{n}{d}}\sum\limits_{g|gcd(d,m)}\varphi(g)\\ &=&\sum\limits_{d|n}\varphi(d)m^{\frac{n}{d}}gcd(d,m) \end{array}$$
对于翻转和颜色平移的 $nm$个变换

1.$n$为奇数，那么翻转轴为一个位置和一个中点构成，此时该位置颜色不能变化，故颜色变换只能为恒同，在翻转下不变的染色方案数为$n\cdot m^{\frac{n+1}{2}}$

2.$n$为偶数，那么翻转轴有两种，第一种是连接两个位置，此时这两个位置颜色也不能变换，颜色变化只能为恒同，在该类变换下不变的染色方案数为$\frac{n}{2}\cdot m^{\frac{n}{2}+1}$；第二种是连接两个中点，此时若$m$为偶数，那么此时颜色变换可以为恒同也可以为颜色平移$\frac{m}{2}$，若$m$为奇数则颜色变换只能恒同，故在该类变换下不变的染色方案数为$gcd(2,m)\cdot \frac{n}{2}\cdot m^{\frac{n}{2}}$

累加以上三个值之后除以$2nm$即可，注意$n$比较大所以要用$Pollard\_rho$算法分解质因数

Code

#include<cstdio>
#include<algorithm>
#include<cmath>
#include<vector>
#include<map>
#include<ctime>
using namespace std;
typedef long long ll;
typedef long double ld;
typedef pair<ll,int>P;
#define maxn 200001
#define mod 998244353
int mul(int x,int y)
{
    ll z=1ll*x*y;
    return z-z/mod*mod;
}
int add(int x,int y)
{
    x+=y;
    if(x>=mod)x-=mod;
    return x; 
}
int Pow(ll x,ll y)
{
    x%=mod,y%=(mod-1);
    int ans=1;
    while(y)
    {
        if(y&1)ans=mul(ans,x);
        x=mul(x,x);
        y>>=1;
    }
    return ans;
}
ll gcd(ll a,ll b)
{
    return b?gcd(b,a%b):a;
}
ll Mul(ll a,ll b,ll c)
{
    return (a*b-(ll)(a/(ld)c*b+1e-3)*c+c)%c;
}
ll Pow(ll a,ll b,ll c)
{
    ll ans=1;
    while(b)
    {
        if(b&1)ans=Mul(ans,a,c);
        a=Mul(a,a,c);
        b>>=1;
    }
    return ans;
}
int cnt=12,prime[]={2,3,5,7,11,13,17,19,23,29,31,37};
bool isprime(ll x)
{
    if(x==2)return 1;
    if(x==1||!(x&1))return 0;
    ll k=x-1,t=0;
    for(;!(k&1);k>>=1)t++;
    for(int i=0;i<cnt&&prime[i]<x;i++) 
    {
        ll g=Pow(prime[i],k,x);
        for(int j=0;j<t;j++)
        {
            ll g0=Mul(g,g,x);
            if(g0==1&&g!=1&&g!=x-1)return 0;
            g=g0;
        }
        if(g!=1)return 0;
    }
    return 1;
}
vector<ll>pf;
vector<P>fact;
ll F(ll t,ll x,ll c)
{
    return (Mul(t,t,x)+c)%x;
}
void pollardRho(ll x)
{
    if(isprime(x)) 
    {
        pf.push_back(x);
        return;
    }
    while(1) 
    {
        ll c=rand()%(x-1)+1,p=rand()%x,q=F(p,x,c);
        while(p!=q) 
        {
            ll d=gcd(x,abs(p-q));
            if(d!=1)
            {
                pollardRho(d);
                pollardRho(x/d);
                return;
            }
            p=F(p,x,c);q=F(F(q,x,c),x,c);
        }
    }
}
void Deal(ll x)
{
    pf.clear();
    fact.clear();
    pollardRho(x);
    sort(pf.begin(),pf.end());
    for(int i=0;i<pf.size();i++)
    {
        int temp=1;
        while(i<pf.size()-1&&pf[i]==pf[i+1])i++,temp++;
        fact.push_back(P(pf[i],temp));
    }
}
int T,ans;
ll n,m;
void dfs(int pos,ll d,ll phi)
{
    if(pos==fact.size())
    {
        ans=add(ans,mul(mul(phi%mod,Pow(m,n/d)),gcd(d,m)%mod));
        return ;
    }
    dfs(pos+1,d,phi);
    ll x=fact[pos].first;
    d*=x,phi*=x-1;
    for(int i=1;i<=fact[pos].second;i++)
    {
        dfs(pos+1,d,phi);
        d*=x,phi*=x;
    }
}
int main()
{
    srand(time(0));
    scanf("%d",&T);
    while(T--)
    {
        scanf("%lld%lld",&n,&m);
        Deal(n);
        ans=0;
        if(n&1)ans=mul(n%mod,Pow(m,n/2+1));
        else 
        {
            ans=mul(n/2%mod,Pow(m,n/2+1));
            ans=add(ans,mul(mul(gcd(m,2),n/2%mod),Pow(m,n/2)));
        }
        dfs(0,1,1);
        ans=mul(ans,Pow(mul(2,mul(n%mod,m%mod)),mod-2));
        printf("%d\n",ans);
    }
    return 0;
}