bzoj 3328 PYXFIB —— 单位根反演

题目：https://www.lydsy.com/JudgeOnline/problem.php?id=3328

单位根反演，主要用到了 \( [k|n] = \frac{1}{k} \sum\limits_{i=0}^{k-1} w_{k}^{in} \)推导见：https://www.cnblogs.com/galaxies/p/bzoj3328.html构造 \( F(x) \) 是为了凑成二项式定理的形式，其实也不难想；又忘记在定义构造函数时写 N 了！把 =N(...) 写成 =(...) 不仅能编译过，还能过样例...以后写构造函数先检查这个！

代码如下：

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
typedef long long ll;
ll rd()
{
  ll ret=0,f=1; char ch=getchar();
  while(ch<'0'||ch>'9'){if(ch=='-')f=0; ch=getchar();}
  while(ch>='0'&&ch<='9')ret=ret*10+ch-'0',ch=getchar();
  return f?ret:-ret;
}
int K,mod,pt[100005],cnt;
ll n;
struct N{
  int a,b,c,d;
  N(int a=0,int b=0,int c=0,int d=0):a(a),b(b),c(c),d(d) {}
  void init(){a=1; b=0; c=0; d=1;}
  N operator * (const N &y) const
  {
    N ret;
    ret.a=((ll)a*y.a+(ll)b*y.c)%mod;
    ret.b=((ll)a*y.b+(ll)b*y.d)%mod;
    ret.c=((ll)c*y.a+(ll)d*y.c)%mod;
    ret.d=((ll)c*y.b+(ll)d*y.d)%mod;
    return ret;
  }
};
ll pw(ll a,ll b){ll ret=1; for(;b;b>>=1,a=a*a%mod)if(b&1)ret=ret*a%mod; return ret;}//ll b
int upt(int x){while(x>=mod)x-=mod; while(x<0)x+=mod; return x;}
int G(int p)
{
  int x=p-1; cnt=0;
  for(int i=2;i*i<=x;i++)
    {
      if(x%i)continue;
      pt[++cnt]=i; if(i*i!=x)pt[++cnt]=x/i;
    }
  for(int a=2;a<p;a++)
    {
      bool fl=1;
      for(int i=1;i<=cnt;i++)
    if(pw(a,pt[i])==1){fl=0; break;}
      if(fl)return a;
    }
}
N pww(N a,ll b)//ll b
{
  N ret; ret.init();
  for(;b;b>>=1,a=a*a)if(b&1)ret=ret*a;
  return ret;
}
int F(int x)
{
  N A=N(x+1,1,1,x); A=pww(A,n);//N()
  x=pw(x,mod-2); x=pw(x,n);
  return (ll)x*A.a%mod;
}
int main()
{
  int T=rd();
  while(T--)
    {
      n=rd(); K=rd(); mod=rd();
      int w=pw(G(mod),(mod-1)/K),inv=pw(w,mod-2);
      int ans=0;
      for(int i=0,t=1;i<K;i++,t=(ll)t*inv%mod)ans=upt(ans+F(t));
      ans=(ll)ans*pw(K,mod-2)%mod;
      printf("%d\n",upt(ans));
    }
  return 0;
}