【2018北京集訓（六）】Lcm

Portal --> 出錯啦qwq（好吧實際上是沒有）

Description

　　給定兩個正整數\(n,k\)，選擇一些互不相同的正整數，知足這些數的最小公倍數剛好爲\(n\)，而且這些數的和爲\(k\)的倍數

　　求選擇的方案數對\(232792561\)取模

　　數據範圍：多組數據，組數\(T<=10,n<=10^{18},k<=20\)，且\(n\)的全部質因子不大於\(100\)

Solution

　　這題。。好神仙啊qwq敲爆腦子都想不出來系列qwq

　　注意到\(n<=10^{18}\)意味着\(n\)至多有\(15\)個不一樣的質因子（前\(15\)個質數乘一下就知道了），而且\(k\)的範圍也是比較小的

　　而後還有一點就是這個模數比較有趣，它\(=lcm(1,2,...,20)+1\)，也就是說它能夠對長度在\(1..20\)的區間DFT

​　　那無論別的咱們能夠先考慮一個最樸素的dp

　　記\(f[i][j]\)表示當前所選的數狀態爲\(i\)，和\(\%k=j\)的方案數，其中這個「所選數的狀態」具體計算方式是：假設寫成分解質因數的形式後\(n=\sum p_i^{mi_i}\)，其中\(p_i\)爲\(n\)的每個不一樣的質因子，咱們的狀態是一個二進制數，而且二進制的每一位對應一個\(n\)的質因子，對於二進制的第\(i\)位（假設對應\(p_i\)），若是說當前所選的數中存在一個數\(x\)知足\(p_i^{mi_i}|x\)，那麼這位爲\(1\)，不然爲\(0\)

​　　這樣一來，咱們最後的答案就應該是\(f[(1<<Cnt)-1][0]\)，其中\(Cnt\)爲\(n\)的不一樣質因子個數

　　至於求解。。分解質因數和預處理因數能夠暴力求解，可是後面的dp直接暴力轉移的話穩穩的T啊（光枚舉子集就要\(3^{Cnt}\)了還要再乘個\(k\))，因此如今的問題是咱們要怎麼比較快速地轉移

​　　觀察這個dp兩維的轉移，咱們會發現第一維的轉移能夠寫成一個集合或卷積的形式（或起來等於某個數就把值累計進去），而第二維的轉移則直接枚舉因數什麼的大力轉移就行了

　　而後這裏有一個頗有意思的想法：首先咱們發現這兩維若是能夠分開處理就會比較好一點，而後注意到這個模數頗有趣，容許咱們進行區間長度\(\in[1,20]\)的DFT，而DFT操做以後，dp的這兩維在某種意義上就是獨立的了，換句話來講，若是說咱們先在同行內轉移\(f\)數組，再對轉之後每行的不徹底計算的\(f\)數組作一次DFT，這行的\(f\)的第二維就不會互相影響了，也就是說咱們能夠經過這種方式實現第一維和第二維轉移的分離

　　這樣咱們就能夠先處理\(f[][i]\)，直觀一點來講就是先在同行轉移\(f\)數組，具體一點就是枚舉\(n\)的因數，而後考慮加進去的貢獻，也就是假設當前枚舉到的因數是\(a[i]\)，\(a[i]\)對應的狀態是\(st\)，那麼咱們能夠對於全部的\(j\in [0,k)\)轉移：
\[ f[st][(j+a[i]\%k)]+=f'[st][j] \]
​　　之因此在後面的\(f\)打了個\('\)是由於咱們這裏要累加進去的是用\(a[i]\)更新前的\(f[st]\)的版本

​　　

​　　這樣咱們就獲得了整合了同行數據以後的\(f\)數組（爲了防止弄混在接下來的描述中咱們仍是將這個不徹底轉移的\(f\)數組記爲\(f1\)好了），而後對於每一行DFT一下，接下來就是考慮第一維的轉移了，更加直觀一點就是。。同列的\(f\)進行轉移

​　　接下來爲了讓描述變得更加簡潔，咱們將第二維省去（由於反正轉移的時候都是同列轉移）

​　　如今咱們要作的事情就是挑出若干個\(f1[i]\)，知足\(i_1\ or\ i_2\ or\ ...\ i_m=st\)而後將這堆\(f1[i]\)的值累加到\(f[st]\)去

　　咱們記\(F(S)=\sum\limits_{i\subseteq S}f[i]\)，若是說咱們知道了\(F(S)\)的取值，那麼只要大力容斥一下就能夠得出\(f\)數組了，具體一點的話就是：
\[ f[T]=\sum\limits_{S\subseteq T}(-1)^{|T|-|S|}F(S) \]
　　咱們將\(F(S)\)進行一下轉化，會發現：
\[ \begin{aligned} F(S)&=\sum\limits_{j\subseteq S}f[j]\\ &=\prod\limits_{i\subseteq S}(1+f1[i]) \end{aligned} \]
　　具體的話就是由於對於\(f\)來講，每個\(f[j]\)都是由若干個\(f1[i]\)組成的，咱們考慮將第二個等號後面的連乘的括號拆掉，展開以後會發現囊括了\(i\)的全部組合方式（爲了方便理解能夠本身用比較小的規模模擬一下），而後由於咱們限制了\(i\subseteq S\)因此能夠保證任意的組合方式都是知足組合後\(\subseteq S\)的

​　　而後因爲\(f1\)是已知的，因此咱們如今就要想如何快速求\(F(S)\)

　　顯然枚舉子集不現實

​　　若是咱們直接從小到大枚舉狀態，每次將這個狀態對應的\(f1\)值轉移的話，可能會出現重複計算的狀況（好比說\(010\)能夠轉移到\(011\)和\(110\)，而\(011\)和\(110\)又均可以轉移到\(111\)，這樣若是按照這種轉移方式的話，\(111\)這裏\(010\)的\(f1\)值就會被重複計算），因此這裏考慮一種很玄學的按順序轉移方式：

​　　咱們考慮二進制一位一位轉移，從低位枚舉到高位，每枚舉一位就把全部的這位爲\(0\)的狀態的值累加到將這位修改成\(1\)後的狀態裏面去，具體一點的話就是：

​　　咱們以\(Cnt=3\)爲例子，狀態轉移大概是這樣：

　　其中橙黃色的箭頭表示轉移的方向，而後位數是從最右邊開始數的

​　　這樣的話顯然全部的狀態都能轉移到它能轉移到的位置，如今咱們來講明一下爲何這樣不會重複計算：首先能夠確認的一點是，這種重複計算的狀況只會出如今轉移到一個被修改了兩個及以上位的狀態的時候，而咱們這樣的枚舉方式每次只轉移到修改了一位的地方，而且是按照數位從低到高，狀態從小到大的順序轉移，好比說剛纔提到的\(010\)在\(111\)中被重複計算的狀況，放在這個轉移方式中就應該是：在轉移第一位的時候，\(010\)會先轉移到\(011\)，而後在轉移第三位的時候再由累加了\(010\)的\(011\)轉移到\(111\)，每次都是從只修改一位的地方轉移過來，不會出現重複計算的狀況（但其實上面的證實也不算。。特別嚴謹。。仍是感性理解既視感qwq）

​　　而後咱們就能夠在\(O(2^{Cnt})\)的時間內求得\(F(S)\)，接下來直接大力容斥一下就能夠得出\(f[(1<<Cnt)-1][0]\)的值最後IDFT一下就能夠得出答案啦

　　

　　代碼大概長這個樣子

#include<iostream>
#include<cstdio>
#include<cstring>
#define ll long long
using namespace std;
const int MOD=232792561,G=71,N=110;
const int p[26]={0,2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97};
ll mi[N],rec_mx[N],Tmp[N],have[N];
ll a[500010];
int f[(1<<15)+10][20];
ll n,K,T,Cnt,all;
void add(int &x,int y){x=(1LL*x+y+MOD)%MOD;}
bool in(int st,int x){return st>>(x-1)&1;}
int St(int x){return 1<<x-1;}
int ksm(int x,int y){
    int ret=1,base=x;
    for (;y;y>>=1,base=1LL*base*base%MOD)
        if (y&1) ret=1LL*ret*base%MOD;
    return ret;
}
namespace DFT{/*{{{*/
    int tmp[110],W[1010][2];
    int inv;
    void init(){
        int rt=ksm(G,(MOD-1)/K);
        W[0][0]=1;
        for (int i=1;i<=1000;++i)
            W[i][0]=1LL*W[i-1][0]*rt%MOD;
        for (int i=0;i<=1000;++i)
            W[i][1]=ksm(W[i][0],MOD-2);
    }
    void dft(int *a,int n,int op){
        for (int i=0;i<n;++i) tmp[i]=0;
        for (int i=0;i<n;++i) 
            for (int j=0;j<n;++j)
                tmp[i]=(1LL*tmp[i]+1LL*a[j]*W[i*j][op==-1]%MOD)%MOD;
        if (op==-1){
            inv=ksm(n,MOD-2);
            for (int i=0;i<n;++i) tmp[i]=1LL*tmp[i]*inv%MOD;
        }
        for (int i=0;i<n;++i) a[i]=tmp[i];
    }
}/*}}}*/
void Div(ll x){
    for (int i=1;i<=25;++i)
        if (x%p[i]==0){
            ++Cnt; rec_mx[Cnt]=1; mi[Cnt]=0; have[Cnt]=p[i];
            while (x%p[i]==0) 
                x/=p[i],++mi[Cnt],rec_mx[Cnt]*=p[i];
        }
}
void dfs(int now,ll prod){
    if (now>Cnt){
        a[++a[0]]=prod;return;
    }
    ll tmp=prod;
    for (int i=0;i<=mi[now];++i){
        dfs(now+1,tmp);
        tmp*=have[now];
    }
}
void prework(){
    a[0]=0; Cnt=0;
    Div(n);
    dfs(1,1);
}
void update(int st,int r){
    for (int i=0;i<K;++i) Tmp[i]=f[st][i];
    for (int i=0;i<K;++i) add(f[st][(i+r)%K],Tmp[i]);
}
void solve(){
    int st,op;
    memset(f,0,sizeof(f));
    all=(1<<Cnt)-1;
    for (int i=0;i<=all;++i) f[i][0]=1;//init
    for (int i=1;i<=a[0];++i){
        st=0;
        for (int j=1;j<=Cnt;++j)
            if (a[i]%rec_mx[j]==0) st|=St(j);
        update(st,a[i]%K);
    }
    for (int i=0;i<=all;++i)
        DFT::dft(f[i],K,1);

    for (int i=1;i<=Cnt;++i)
        for (int j=0;j<=all;++j)
            if (in(j,i))
                for (int k=0;k<K;++k)
                    f[j][k]=1LL*f[j][k]*f[j^St(i)][k]%MOD;
    for (int i=0;i<all;++i){
        op=1;
        for (int j=1;j<=Cnt;++j)
            if (!in(i,j)) op*=-1;
        for (int j=0;j<K;++j)
            add(f[all][j],op*f[i][j]);
    }
    DFT::dft(f[all],K,-1);
    printf("%d\n",f[all][0]);
}

int main(){
#ifndef ONLINE_JUDGE
    freopen("a.in","r",stdin);
#endif
    scanf("%d",&T);
    for (int o=1;o<=T;++o){
        scanf("%lld%lld",&n,&K);
        DFT::init();
        prework();
        solve();
    }
}