单位根反演

　单位根反演

感觉这东西挺冷门的......

单位根反演

总之先记个柿子，对于整数\(n,k\)，有
\[ [n \mid k] = \frac{1}{n}\sum\limits_{i = 0} ^ {n-1} \omega_{n}^{ik} \]
证明：

• 如果\(n \mid k\)，那么\(\omega_{n}^{ik} = \omega_{n}^0 =1\)，原式乘上\(\frac{1}{n}\)后显然是\(1\)。

• 如果\(n \not\mid k\)，那么等比数列求和
  \[ \sum\limits_{i = 0} ^ {n-1} \omega_n^{ik} = \sum\limits_{i = 0} ^ {n-1} (\omega_{n}^{k})^i = \frac{\omega_{n}^{kn} - \omega_{n}^{0}}{\omega_{n}^{k} - 1} = 0 \]

没了。

例题

bzoj3328 PYXFIB

题面太妙了......这里在说一遍题意

给整数\(n,k,p\)，让你求
\[ \sum\limits_{i=0}^{\left\lfloor \frac{n}{k} \right\rfloor} \binom{n}{ik} f_{ik} \ mod \ p \]
其中，\(1\le n\le 10^{18}, 1\le k \le 20000\)，保证\(p\)是质数且\(p \equiv 1 \ (mod \ k)\)，\(f_i\)表示斐波那契数列的第\(i\)项（\(f_0 = 1, f_1 = 1, f_i = f_{i-1}+f_{i-2}(i>1)\)）。

我们先把下取整去掉
\[ \sum\limits_{i = 0}^n [k \mid i] \binom{n}{i}f_i \]
组合数很想直接二项式定理......但我们有几个问题

• 1. 后面是斐波那契数列
• 1. 还有一个整除的限制

斐波那契好解决，直接用矩乘，即\(A = \begin{bmatrix} 1 & 1\\ 1 & 0 \end{bmatrix}\)，那么\(f_n = (A^n)_{1,2} + (A^n)_{2,2}\)，即\(A^n\)的最后一列上两个数的和。

下面不看第\(2\)个限制，就是一个裸的二项式定理（二项式定理对矩阵也使用哦）
\[ \sum\limits_{i = 0} ^ n \binom{n}{i} f_i = ((A + I)^n)_{1,2} + ((A+I)^n)_{2,2} \]
那前面整除咋办呢？注意到\(k\)只有\(20000\)，单位根反演暴力拆
\[ \begin{aligned}&\quad \sum\limits_{i = 0}^n \frac{1}{k}\sum\limits_{j=0}^{k} \omega_{k}^{ij} \binom{n}{i}f_i \\&= \frac{1}{k} \sum_{j = 0} ^ k \sum\limits_{i = 0} ^ n \binom{n}{i} (\omega_{k}^{j})^if_i\end{aligned} \]
然后直接把\(\omega_{k}^{j}\)乘到矩阵\(A\)里，再二项式定理
\[ \frac{1}{k}\sum\limits_{j=0}^k ((A\omega_{k}^j + I)^n)_{1,2} + ((A\omega_{k}^j + I)^n)_{2,2} \]
对\(p\)取模的话直接用原根代替单位根，再做矩阵快速幂就好了，复杂度\(O(k \ log \ n)\)。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
int mod,K,g;ll n;
inline int fpow(ll a,ll b){
    int ret=1; for (a%=mod;b;b>>=1,a=1ll*a*a%mod)
        if (b&1) ret=1ll*ret*a%mod;
    return ret;
}
inline int check(int k1){
    for (int i=2;i*i<=mod-1;i++)
        if ((mod-1)%i==0&&(fpow(k1,i)==1||fpow(k1,(mod-1)/i)==1))return 0;
    return 1;
}
inline void getG(){
    for (int i=2;i<=100;i++)
        if (check(i)){g=i;break;}
}
inline int add(int x,int y){return x+y>=mod?x+y-mod:x+y;}
struct matrix{
    int a[3][3];
    matrix(int R=2,int C=2){
        for (int i=0;i<R;i++)
            for (int j=0;j<C;j++) a[i][j]=0;
    }
    int* operator [] (int k1){return a[k1];}
    matrix operator * (matrix k1)const{
        matrix res;
        for (int i=0;i<2;i++)
            for (int j=0;j<2;j++)
                for (int k=0;k<2;k++)res[i][j]=add(res[i][j],1ll*a[i][k]*k1[k][j]%mod);
        return res;
    }
    matrix operator * (ll k1)const{
        k1%=mod; matrix ret;
        for (int i=0;i<2;i++)
            for (int j=0;j<2;j++)
                ret[i][j]=1ll*a[i][j]*k1%mod;
        return ret;
    }
    matrix operator + (matrix k1)const{
        matrix ret;
        for (int i=0;i<2;i++)
            for (int j=0;j<2;j++) ret[i][j]=add(a[i][j],k1[i][j]);
        return ret;
    }
    void print(){
        for (int i=0;i<2;i++,puts(""))
            for (int j=0;j<2;j++)
                printf("%d ",a[i][j]);
    }
};
inline matrix fmpow(matrix A,ll t){
    matrix ret; ret[0][0]=ret[1][1]=1;
    for (;t;t>>=1,A=A*A) if (t&1) ret=ret*A;
    return ret;
}
int main(){
    matrix I; I[0][0]=I[1][1]=1;
    matrix A; A[0][0]=A[0][1]=A[1][0]=1;
    int T;scanf("%d",&T);while(T--){
        scanf("%lld%d%d",&n,&K,&mod);
        getG(); int ans=0,wk=fpow(g,(mod-1)/K);
        for (int i=0,w=1;i<K;i++,w=1ll*w*wk%mod){
            matrix B=fmpow(A*w+I,n); //B.print();
            ans=add(ans,B[0][1]),ans=add(ans,B[1][1]);
        }
        ans=1ll*ans*fpow(K,mod-2)%mod;
        printf("%d\n",ans);
    }
    return 0;
}

还有题目先鸽着吧......看心情会更的......

单位根反演

原文：https://www.cnblogs.com/wxq1229/p/12242177.html