基础数论刷题柱
题单来自于：[2025年暑假集训B班基础数论作业]
jzyz1611 Interval GCD
如果用线段树维护区间的gcd，难以修改
注意到gcd(a,b)=gcd(a,a-b)，因此区间[l,r]的gcd可以修改为询问gcd(a[l],a[l]-a[l+1],a[l+1]-a[l+2]...a[r-1]-a[r])，此时询问a[l]是区间加，单点询问，非常简单
设b[i]是a[i]的差分数组，也就是b[i]=a[i]-a[i-1]，维护b[i]就可以单点修改区间询问，也不困难
对于区间加，单点询问的线段树可以不下传懒标记，而是把区间加放在区间节点上，单点询问时把路径上的区间加操作加在一起。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=500010,mod=998244353;//1e9+7;
ll read(){ll x;scanf("%lld",&x);return x;}
ll quick(ll a,ll b){ll t=1;while(b){if(b&1)t=t*a%mod;b=b/2;a=a*a%mod;}return t;}
int lowbit(int x){return x&(-x);}
ll gcd1(ll a,ll b)
{
a=abs(a);
b=abs(b);
if(b==0)
return a;
return gcd1(b,a%b);
}
int n,q;
ll a[N],b[N],sum[N*4],lz[N*4],GCD[N*4];
void build(int x,int l,int r)
{
if(l==r)
{
GCD[x]=b[l];
sum[x]=a[l];
return ;
}
int mid=(l+r)/2;
build(x*2,l,mid);
build(x*2+1,mid+1,r);
GCD[x]=gcd1(GCD[x*2],GCD[x*2+1]);
}
void dan(int x,int l,int r,int pos,ll v)
{
if(l==r)
{
GCD[x]+=v;
return ;
}
int mid=(l+r)/2;
if(pos<=mid)
dan(x*2,l,mid,pos,v);
else
dan(x*2+1,mid+1,r,pos,v);
GCD[x]=gcd1(GCD[x*2],GCD[x*2+1]);
}
void qu(int x,int l,int r,int tl,int tr,ll v)
{
if(tl<=l&&r<=tr)
{
lz[x]+=v;
return ;
}
int mid=(l+r)/2;
if(tl<=mid)
qu(x*2,l,mid,tl,tr,v);
if(tr>mid)
qu(x*2+1,mid+1,r,tl,tr,v);
}
ll askqu(int x,int l,int r,int tl,int tr)
{
if(tl<=l&&r<=tr)
return GCD[x];
int mid=(l+r)/2;
ll t=0;
if(tl<=mid)
t=gcd1(t,askqu(x*2,l,mid,tl,tr));
if(tr>mid)
t=gcd1(t,askqu(x*2+1,mid+1,r,tl,tr));
return t;
}
//有的时候不需要下传
//把修改挂在完全覆盖的，最少的区间上
//询问时，路径上的lz相加，再加上叶子的
ll askdan(int x,int l,int r,int pos)
{
if(l==r)
return sum[x]+lz[x];
int mid=(l+r)/2;
if(pos<=mid)
return askdan(x*2,l,mid,pos)+lz[x];
else
return askdan(x*2+1,mid+1,r,pos)+lz[x];
}
int main()
{
// freopen("1.in","r",stdin);
n=read();q=read();
for(int i=1;i<=n;i++)
a[i]=read();
for(int i=1;i<n;i++)
b[i]=a[i]-a[i+1];
build(1,1,n);
for(;q;q--)
{
char c;
cin>>c;
if(c=='C')
{
int l=read(),r=read();ll x=read();
if(l!=1)
dan(1,1,n,l-1,-x);
dan(1,1,n,r,x);
qu(1,1,n,l,r,x);
}
else
{
int l=read(),r=read();
if(l==r)
cout<<askdan(1,1,n,l)<<'\n';
else
cout<<abs(gcd1(askqu(1,1,n,l,r-1),askdan(1,1,n,l)))<<'\n';
}
}
}
//1、以1为根节点
//2、x的左儿子是 2x，右儿子是 2x+1
//3、x是[l,r]，如果l==r叶子没有儿子
//mid=(l+r)/2
//左儿子:l,mid
//右儿子:mid+1,r
//4、线段树4*n空间
//5、区间修改懒标记
//sum[x]+=(r-l+1)*v
//lz[x]+=v
//表示儿子们应该+v但是还没加
jzyz1396 [CQOI]余数之和-整除分块
\(k\mod i=k-\lfloor k/i\rfloor *i\)，因此可以把\(\sum k \mod i\)变为\(n*k- \sum \lfloor k/i\rfloor*i\)
在整除分块时，lr区间内\(k/i\)的值相同，因此可以加速
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll n,k,ans;
int main()
{
cin>>n>>k;
ans=n*k;
for(ll l=1,r;l<=n&&l<=k;l=r+1)
{
r=k/(k/l);
r=min(r,n);
r=min(r,k);
ans=ans-(k/l)*(l+r)*(r-l+1)/2;
}
cout<<ans;
}
河南萌新联赛2025第（二）场：河南农业大学 A 约数个数和
\(\sum d(i)\)可以转化为\(\sum n/i\)，于是变成了整除分块模板
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll n,k,ans;
int main()
{
cin>>n>>k;
ans=n*k;
for(ll l=1,r;l<=n&&l<=k;l=r+1)
{
r=k/(k/l);
r=min(r,n);
r=min(r,k);
ans=ans-(k/l)*(l+r)*(r-l+1)/2;
}
cout<<ans;
}
BZOJ4488 [Jsoi2015]最大公约数
如果固定右左端点，则区间gcd随着右端点向右移动而减小，每次减小必然除以二或者除的更多，因此突变点最多有log个
接下来只需移动左端点时更新突变点和区间gcd，枚举和gcd更新答案即可。维护突变点稍有难度，考虑左端点从大往小变，这样原来不是突变点的依然不是突变点，因为原来不是突变点说明\(\gcd(a[l+1],a[l+2]...a[pos-2],a[pos-1])|a[pos]\)，那么加上a[l]后，依然满足\(\gcd(a[l],a[l+1],a[l+2]...a[pos-2],a[pos-1])|a[pos]\)，依然不会是突变点。如果原来是突变点，那么加上\(a[l]\)后可能是突变点也可能不是。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=500010,mod=998244353;//1e9+7;
ll read(){ll x;scanf("%lld",&x);return x;}
ll quick(ll a,ll b){ll t=1;while(b){if(b&1)t=t*a%mod;b=b/2;a=a*a%mod;}return t;}
int lowbit(int x){return x&(-x);}
int n;
ll a[N],ans;
vector<pair<int,ll>>o[N];
int main()
{
n=read();
for(int i=1;i<=n;i++)
a[i]=read();
for(int l=n;l>=1;l--)
{
o[l].push_back({l,a[l]});
for(auto [pos,v]:o[l+1])
{
if(__gcd(v,a[l])==o[l].back().second)
continue;
else
o[l].push_back({pos,__gcd(v,a[l])});
}
for(int j=0;j<o[l].size();j++)
{
if(j==o[l].size()-1)
ans=max(ans,o[l][j].second*(n-l+1));
else
ans=max(ans,o[l][j].second*(o[l][j+1].first-l));
}
}
cout<<ans;
}
//随着r增大
//gcd要么不变，要么变小
//变小至少除以二
//最多变小log(1e12)次
//二分突变点不太行
//已知l+1的全体突变点
//怎么拿到l的全体突变点
//pos=l
//原来突变的，这次有可能突变
//a[pos]不整除gcd(l+1,pos-1)
//a[pos] gcd(l,pos-1)
//原来不突变的，这次一定不突变
//gcd(l+1,pos-1)|a[pos]
//gcd(l,pos-1)|gcd(l+1,pos-1)
//--->
//gcd(l,pos-1)|a[pos]
p2492. 习题10.2.4 质因数分解
\(n\)只有\(10^{12}\)，使用暴力的\(t*\sqrt{10^{12}}\)即可。
如果预处理质数，用质数分解质因数会更快一点。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=500010,mod=998244353;//1e9+7;
ll read(){ll x;scanf("%lld",&x);return x;}
ll quick(ll a,ll b){ll t=1;while(b){if(b&1)t=t*a%mod;b=b/2;a=a*a%mod;}return t;}
int lowbit(int x){return x&(-x);}
void work()
{
ll n=read();
for(ll i=2;i<=n/i;i++)
{
while(n%i==0)
{
n=n/i;
cout<<i<<' ';
}
}
if(n!=1)
cout<<n<<' ';
cout<<'\n';
}
int main()
{
freopen("prime.in","r",stdin);
freopen("prime.out","w",stdout);
for(int t=read();t;t--)
work();
}
loj2314. 「NOIP2017」小凯的疑惑
证明我也不会，结论是输出\(a*b-a-b\)。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=1000010,mod=998244353;//1e9+7;
ll read(){ll x;scanf("%lld",&x);return x;}
ll quick(ll a,ll b){ll t=1;while(b){if(b&1)t=t*a%mod;b=b/2;a=a*a%mod;}return t;}
int lowbit(int x){return x&(-x);}
int main()
{
ll a=read(),b=read();
cout<<a*b-a-b;
}
jzyz381. 【bzoj2190】[SDOI2008]仪仗队
1、可以发现\((2,6)\)不能被看见，因为被\((1,3)\)挡住了，可以发现一般情况下需要\(gcd(x,y)\neq 1\)才能被看见，否则\((x,y)\)会被\((x/gcd,y/gcd)\)挡住。
特殊情况：除了\((0,0)\)以外不存在\(x=y\)的了，除了\((0,1),(1,0)\)以外不存在\(x=0\)或\(y=0\)的。
2、\(3\)个特殊的拿走后，做一下对称，答案是满足\((x\geq 1,y>1,gcd(x,y)=1,x<y)\)的数量\(*2\)再\(+3\)，这就是\(\varphi(2)+\varphi(3)+...\varphi(n)\)
3、\(n\)很小，暴力求欧拉函数也可以
\(1\sim n\)把包含\(p1\)的删掉
\(n*(1-1/p1)\)
再把含p2的删掉
\(n*(1-1/p1)*(1-1/p2)\)
再把含\(p3\)的删掉
\(n*(1-1/p1)*(1-1/p2)*(1-1/p3)\)
...
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=1000010,mod=998244353;//1e9+7;
ll read(){ll x;scanf("%lld",&x);return x;}
ll quick(ll a,ll b){ll t=1;while(b){if(b&1)t=t*a%mod;b=b/2;a=a*a%mod;}return t;}
int lowbit(int x){return x&(-x);}
int n,ans;
int phi(int x)
{
int res=x;
for(int i=2;i<=x/i;i++)
if(x%i==0)
{
res=res/i*(i-1);
while(x%i==0)
x=x/i;
}
if(x!=1)
res=res/x*(x-1);
return res;
}
int main()
{
freopen("queue.in","r",stdin);
freopen("queue.out","w",stdout);
n=read()-1;
for(int i=2;i<=n;i++)
ans=ans+phi(i);
ans=ans*2+3;
cout<<ans;
}
BZOJ2705 洛谷P2303 Longge的问题
\(gcd(i,n)\)，\(gcd\)结果一定是\(n\)的因子。
如果固定\(gcd\)为\(d\)，\(d\)的贡献是\(d*\)有几个数字\(i\)和\(n\)的\(gcd\)恰好为\(d\)。
拿着\(gcd(i,n)=d,1\leq i \leq n\)开始思考：同时除以\(d\)，问题转变成\(gcd(i,n/d)=1,1\leq i \leq n/d\)，这其实就是\(\varphi (n/d)\)。
结论是暴力枚举因子\(d\)，求\(\sum d*\varphi(n/d)\)。
注意n需要开ll
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=1000010,mod=998244353;//1e9+7;
ll read(){ll x;scanf("%lld",&x);return x;}
ll quick(ll a,ll b){ll t=1;while(b){if(b&1)t=t*a%mod;b=b/2;a=a*a%mod;}return t;}
int lowbit(int x){return x&(-x);}
ll n,ans;
ll phi(ll x)
{
ll res=x;
for(int i=2;i<=x/i;i++)
if(x%i==0)
{
res=res/i*(i-1);
while(x%i==0)
x=x/i;
}
if(x!=1)
res=res/x*(x-1);
return res;
}
int main()
{
n=read();
for(int d=1;d<=n/d;d++)
{
if(n%d)
continue;
ans=ans+1ll*d*phi(n/d);
if(n/d!=d){
ans=ans+1ll*n/d*phi(d);
}
}
cout<<ans;
}
BZOJ3884. 上帝与集合的正确用法
设 \(ans=2^{2^{2^{2...}}}\)
则根据扩展欧拉定理，\(ans \mod p=2^{ans\mod \varphi(p)+\varphi (p)}\mod p\)
左边是要求的，右边只有\(ans \mod \varphi(p)\)不知道。递归地去求\(ans\mod \varphi(p)\)，依然使用拓欧，直到模数为\(1\)，此时可以直接返回\(0\)。
在oj上跑了1.35s，在洛谷上跑了300ms，大胆推测oj是洛谷的四倍（）
使用记忆化可以使得平均下降约200ms，但是我在oj上调大了时限，都能过。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=1000010;//1e9+7;
ll read(){ll x;scanf("%lld",&x);return x;}
ll quick(ll a,ll b,ll mod){ll t=1;while(b){if(b&1)t=t*a%mod;b=b/2;a=a*a%mod;}return t;}
int lowbit(int x){return x&(-x);}
int cnt,pri[700000],v[10000010],phi[10000010];
int ask(int p)
{
if(p==1)
return 0;
return quick(2,ask(phi[p])+phi[p],p);
}
void work()
{
cout<<ask(read())<<'\n';
}
int main()
{
// freopen("queue.in","r",stdin);
// freopen("queue.out","w",stdout);
phi[1]=1;
for(int i=2;i<=10000000;i++)
{
if(v[i]==0){
pri[++cnt]=i;
v[i]=i;
phi[i]=i-1;
}
for(int j=1;1ll*i*pri[j]<=10000000&&j<=cnt&&pri[j]<=v[i];j++)
{
v[i*pri[j]]=i;
if(i%pri[j]==0)//有多个pri[j]
phi[i*pri[j]]=phi[i]*pri[j];
else
phi[i*pri[j]]=phi[i]*(pri[j]-1);
}
}
for(int t=read();t;t--)
work();
}
BZOJ2186 [Sdoi2008] 沙拉公主的困惑
求\(1\sim n!\)里和\(m!\)互质的数量。
由于\(gcd(i,m!)=gcd(i-m!,m!)\)，因此把\(1\sim n!\)分成\(n!/m!\)份，每一份的互质情况都一样，求\(1\sim m!\)里有几个和\(m!\)互质，最后乘一下就好了。
\(1\sim m!\)里有几个和\(m!\)互质，这就是欧拉函数的定义，根据\(\varphi(n)=n\prod{(p_i-1)\over p_i}\)，这里的\(p_i\)包含了\(1\sim m\)里的所有质数。由于多组询问，可以把这个\(\prod\)前缀预处理一下。
20250812更新：类似于组合数的大\(n,m\)小模数，本题如果模数过小也会遭遇一些问题，比如\(n!/m!*phi(m!)\)虽然\(\%mod\)不是\(0\)，但是\(n!\)和\(m!\)和\(\varphi(m!)\)可能有多个\(mod\)，为\(0\)了。我们掏出拓展卢卡斯的做法，对\(n!\)和\(\varphi(i)\)都处理一下包含了几个\(mod\)，把\(mod\)的数量统计一下，和\(mod\)互质的部分乘在一起。这样输出的时候判断一下\(mod\)的数量是不是\(0\)就可以知道答案是不是\(0\)。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=10000010;
int phi[N],f[N],mod=998244353;
ll quick(ll a,ll b){ll t=1;while(b){if(b&1)t=t*a%mod;a=a*a%mod;b=b/2;}return t;}
int lowbit(int x){return x&(-x);}
ll read(){ll x;scanf("%lld",&x);return x;}
int pri[N],v[N],cnt,numphi[N],numf[N];
int main()
{
// freopen("1.in","r",stdin);
// freopen(".out","w",stdout);
int t=read();
mod=read();
for(int i=2;i<=10000000;i++)
{
if(v[i]==0){
pri[++cnt]=i;
v[i]=i;
}
for(int j=1;1ll*i*pri[j]<=10000000&&j<=cnt&&pri[j]<=v[i];j++)
v[i*pri[j]]=i;
}
phi[0]=phi[1]=f[0]=f[1]=1;
for(int i=2;i<=10000000;i++)
{
if(v[i]==i){
int t=(i-1);
while(t%mod==0)t=t/mod,numphi[i]++;
phi[i]=1ll*phi[i-1]*t%mod;
numphi[i]+=numphi[i-1];//numphi[i]表示phi[i]
//里面mod的数量
}
else{
int t=i;
while(t%mod==0)t=t/mod,numphi[i]++;
phi[i]=1ll*phi[i-1]*t%mod;//phi[i!]
numphi[i]+=numphi[i-1];
}
int t=i;
while(t%mod==0)t=t/mod,numf[i]++;
f[i]=1ll*f[i-1]*t%mod;//预处理阶乘
numf[i]+=numf[i-1];
}
for(;t;t--)
{
int n=read(),m=read();
if(numf[n]-numf[m]+numphi[m]>0)
cout<<"0\n";
else
cout<<f[n]*quick(f[m],mod-2)%mod*phi[m]%mod<<'\n';
}
}
BZOJ2749 [HAOI2012] 洛谷P2350 外星人
给你一个\(N\)，由于\(N\)很大，以唯一分解的形式给你。
设\(\varphi^x(N)\)表示\(\varphi(\varphi(\varphi(...\varphi(N))\)，迭代\(x\)次。
手玩几个可以发现2很重要
如果刚开始是\(3^3\)
那么下一次\(2^1*3^2\)
下一次\(3^1*2^1\)
下一次\(2^1\)
还需要一次
如果刚开始是\(5^3\)
下一次\(2^2*5^2\)
下一次\(2^3*5^1\)
下一次\(2^4\)
还需要四次
一个\(5\)=两个\(2\)，因此答案是\(2*3+1\)
如果3和5同时来呢
答案是\(1+q[3]+q[5]*2\)其中\(q[i]\)表示\(i\)的幂次而不是题意的\(q[i]\)
如果2和3和5同时来呢
答案是\(q[2]+q[3]+q[5]*2\)
这个系数到底是啥呢，比如\(11\)会变成\(10=2*5\)，\(2\)和\(5\)分开计算，分别会贡献\(1\)个\(2\)和\(2\)个\(2\)。设\(f[i]\)表示\(i\)对\(2\)的贡献，如果i是质数，\(f[i]=f[i-1]\)，否则\(f[i]=\sum f[p_j]*k_j\)其中\(p_j\)和\(k_j\)是i的质因数分解后的结果。可以发现你任选一个\(p_j\)，\(f[i]\)也等于\(f[i/p_j]+f[p_j]\)这个可以线性筛的时候计算一下。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=10000010;
ll mod=998244353;
ll quick(ll a,ll b){ll t=1;while(b){if(b&1)t=t*a%mod;a=a*a%mod;b=b/2;}return t;}
int lowbit(int x){return x&(-x);}
ll read(){ll x;scanf("%lld",&x);return x;}
int pri[N],v[N],f[N],cnt;
int main()
{
// freopen(".in","r",stdin);
// freopen(".out","w",stdout);
f[1]=1;
for(int i=2;i<=100000;i++)
{
if(v[i]==0){
pri[++cnt]=i;
v[i]=i;
f[i]=f[i-1];
}
for(int j=1;1ll*i*pri[j]<=100000&&j<=cnt&&pri[j]<=v[i];j++){
v[i*pri[j]]=i;
f[i*pri[j]]=f[i]+f[pri[j]];
}
}
for(int t=read();t;t--)
{
ll ans=1;//假装没有2
for(int n=read();n;n--)
{
int p=read();
if(p==2)//发现有2，此时让ans--
ans--;
ans=ans+f[p]*read();
}
cout<<ans<<'\n';
}
}
BZOJ1876 [SDOI2009] SuperGCD
高精度，求gcd
python秒了（bushi）
import math
import sys
sys.set_int_max_str_digits(10010)
a=int(input())
b=int(input())
print(math.gcd(a,b))
使用辗转相除法很不妙，一个是容易超时，一个是高精度除法不好写。
使用普通更相减损术，复杂度可能是\(O(n)\)的，比如一个大奇数和一个小\(2\)，一直减一直减...
考虑奇偶加速的更相减损术：如果\(ab\)都是偶数，可以同时除以二，令计数器\(cnt++\)。如果其中一个是偶数，那么可以把他除以二。如果两个都是奇数，则可以大的\(-=\)小的。
普通高精度复杂度吃满了\(O(n*log_2(10^n))\approx 3n^2\)，本地跑第四个测试点需要\(3.5s\)，我写了个压位高精度才能通过。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=10000010;
ll mod=998244353;
ll quick(ll a,ll b){ll t=1;while(b){if(b&1)t=t*a%mod;a=a*a%mod;b=b/2;}return t;}
int lowbit(int x){return x&(-x);}
ll read(){ll x;scanf("%lld",&x);return x;}
void input(vector<int>&a)
{
string s;
cin>>s;
reverse(s.begin(),s.end());
for(int i=0;i*4<s.size();i++)
{
a.push_back(0);
for(int j=min(i*4+3,(int)s.size()-1);j>=i*4;j--)
a.back()=a.back()*10+s[j]-'0';
}
}
bool operator <(vector<int>&a,vector<int>&b)
{
if(a.size()<b.size())
return 1;
if(a.size()>b.size())
return 0;
for(int i=a.size()-1;i>=0;i--)
{
if(a[i]<b[i])
return 1;
if(a[i]>b[i])
return 0;
}
return 0;//相等return 0
}
void operator -=(vector<int>&a,vector<int>&b)
{
for(int i=0;i<b.size();i++)
a[i]-=b[i];
for(int i=0;i<a.size();i++)
if(a[i]<0){
a[i+1]--;
a[i]+=10000;
}
while(a.size()&&a.back()==0)
a.pop_back();
}
void operator /=(vector<int>&a,int v)
{
for(int i=a.size()-1;i>=0;i--)
{
if(a[i]%2==0)
a[i]=a[i]/2;
else
a[i]=a[i]/2,a[i-1]+=10000;
}
if(a.back()==0)
a.pop_back();
}
void operator *=(vector<int>&a,int v)
{
a.push_back(0);
for(int i=0;i<a.size();i++)
a[i]=a[i]*2;
for(int i=0;i<a.size();i++)
if(a[i]>=10000)
{
a[i]-=10000;
a[i+1]++;
}
if(a.back()==0)
a.pop_back();
}
int main()
{
// freopen(".in","r",stdin);
// freopen(".out","w",stdout);
vector<int>a,b;
input(a);
input(b);
int cnt=0,num=0;
while(a.size()!=0&&b.size()!=0)
{
if(a<b)
swap(a,b);
num++;
if(a[0]%2==0&&b[0]%2==0)
{
cnt++;
a/=2;
b/=2;
}
else if(a[0]%2==0){
a/=2;
}
else if(b[0]%2==0)
b/=2;
else
a-=b;
}
if(a.size()==0)
a=b;
for(int i=1;i<=cnt;i++)
a*=2;
reverse(a.begin(),a.end());
cout<<a[0];
for(int i=1;i<a.size();i++)
printf("%04d",a[i]);
}
洛谷P2480 古代猪文
样例是\(2^{C(4,1)+C(4,2)+C(4,4)}=2^{11}=2048\)
模数\(999911659\)是个质数，根据扩展欧拉定理，我们求$$sum=\sum_{d|n} C(n,d)\ mod\ 999911658 $$，再输出\(p^{sum} \% 999911659\)即可
\(999911658=2*3*4679*35617\)，质因数幂次都是\(1\)，可以不用拓展卢卡斯定理来做。可以对于每个\(C(n,d)\)手动求出对四个质数的余数，用中国剩余定理合并起来。
以模\(4679\)意义下的组合数为例，是不能用阶乘和阶乘的逆元来求的，你需要使用卢卡斯定理处理之。
需要特判\(g=mod\)的情况。因为如果恰好\(mod\)也是\(0\)的话，根据我们的拓展欧拉定理，应该输出\(p^{sum+\varphi(999911659)}=0\)，不写的话是\(p^0\)为1。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll quick(ll a,ll b,int mod=999911659){ll t=1;while(b){if(b&1)t=t*a%mod;a=a*a%mod;b=b/2;}return t;}
const int N=10000010,mod=999911659;
int lowbit(int x){return x&(-x);}
ll read(){ll x;scanf("%lld",&x);return x;}
int m[5],lx[5][36000],inv[5][36000],a[5],M[5],t[5],n,g,sum;
int lucas(int a,int b,int i)
{
if(a<b)
return 0;
if(a<m[i]&&b<m[i])
return lx[i][a]*inv[i][b]%m[i]*inv[i][a-b]%m[i];
return lucas(a%m[i],b%m[i],i)*lucas(a/m[i],b/m[i],i)%m[i];
}
int ask(int k)
{
int x=0;
for(int i=1;i<=4;i++){
a[i]=lucas(n,k,i);
// printf("C(%d,%d)%%%d=%d\n",n,k,m[i],lucas(n,k,i));
x=(x+1ll*a[i]*M[i]%(mod-1)*t[i])%(mod-1);
}
return x;
}
int main()
{
// freopen(".in","r",stdin);
// freopen(".out","w",stdout);
m[1]=2;
m[2]=3;
m[3]=4679;
m[4]=35617;
for(int i=1;i<=4;i++)
{
M[i]=(mod-1)/m[i];
for(t[i]=1;t[i]*M[i]%m[i]!=1;t[i]++)
;
}
for(int i=1;i<=4;i++)
{
lx[i][0]=inv[i][0]=1;
for(int j=1;j<m[i];j++){
lx[i][j]=lx[i][j-1]*j%m[i];
inv[i][j]=quick(lx[i][j],m[i]-2,m[i]);
}
}
n=read(),g=read();
for(int k=1;k<=n/k;k++)
{
if(n%k)
continue;
sum=(sum+ask(k))%(mod-1);
if(n/k!=k)
sum=(sum+ask(n/k))%(mod-1);
}
if(g==mod)
cout<<0;
else
cout<<quick(g,sum);
}
我们C班数论进阶再见！ヾ(￣▽￣)ByeBye
浙公网安备 33010602011771号