同余系基本定理2

command_block · 2020-04-04 12:57:04 · 个人记录

上回请看同余系基本定理1

默认大家已经对同余系有一定的了解,而且能够熟练运用基础算法。

P1495 【模板】中国剩余定理(CRT)

有若干个同余方程\begin{cases}x=c_1\pmod {p_1}\\x=c_2\pmod {p_2}\\...\\x=c_m\pmod {p_m}\end{cases}

满足p_{1...m}两两互质,求x的最小正整数解。

考虑把这些方程两两合并。

先看两个方程\begin{cases}x=c_1\pmod{p_1}\\x=c_2\pmod {p_2}\end{cases}

可以构造x=c_1p_2+c_2p_1,这样在模p_1时显出c_1p_2,模p_2时显出c_2p_1.

问题在于,我们希望得到c_1而非c_1p_2。这好办,再乘上p_2在模p_1时的逆即可。

合并完之后的模数是p_1p_2,由于互质所以总能求逆。

边界方程可以视作x=0\pmod 1.

复杂度O(m\log p).

#include<cstdio>
#define ll long long
using namespace std;
void exgcd(ll a,ll b,ll &x,ll &y){
  if (b==0){x=1;y=0;return ;}
  exgcd(b,a%b,y,x);y-=(a/b)*x;
}
ll inv(ll a,ll m){
  ll x,y;exgcd(a,m,x,y);
  return (x%m+m)%m;
}
int n;
ll c,p,ret,mp;
int main()
{
  scanf("%d",&n);
  ret=0;mp=1;
  for (int i=1;i<=n;i++){
    scanf("%lld%lld",&p,&c);
    ll sp=mp*p;
    ret=(ret*inv(p,mp)%sp*p+mp*inv(mp,p)%sp*c)%sp;
    mp=sp;
  }printf("%lld",ret);
  return 0;
}

例题 : P3868 [TJOI2009]猜数字
P3807 【模板】卢卡斯定理

$$\dbinom{n}{m}\bmod p=\dbinom{\lfloor n/p\rfloor}{\lfloor m/p\rfloor}\dbinom{n\bmod p}{m\bmod p}\bmod p$$ $\color{blue}\text{证明:}$ (可不掌握) 注意到$[x^m](1+x)^n=\dbinom{n}{m}

构造(1+x)^p=1+\binom{p}{1}x+\binom{p}{2}x+...+\binom{p}{p}x^p

对于1\leq i\leq p-1,总有\binom{p}{i}=0\pmod p.

原因是\dbinom{p}{i}=\dfrac{p!}{(p-i)!i!},而p是素数,所以 (p-i)! 和 i! 都不含因子p,无法抵消p!中的恰一个因子p.

那么可以得到(1+x)^p=1+x^p\pmod p

能进一步推出(1+x)^{p^k}=(1+x^p)^{p^{k-1}}=...=1+x^{p^k}

我们把n,m都表示成两个p进制数,即\begin{cases}n=\sum\limits_{i=0}a_ip^i\\m=\sum\limits_{i=0}b_ip^i\end{cases}

然后考虑(1+x)^n=\prod\limits_{i=0}(1+x)^{a_ip^i}

(1+x)^n=\prod\limits_{i=0}(1+x^{p^k})^{a^i}\pmod p

提取第m项系数。

[x^m](1+x)^n=[x^m]\prod\limits_{i=0}(1+x^{p^k})^{a^i}\pmod p

由于乘积项的次数不交,可以把m分解。

“不交”的意思是某两位之间不会互相影响。因为将一个 p 进制数的后 n 位加起来之后肯定不足第 n+1 位的一个 1 大。

\dbinom{n}{m}=\prod\limits_{i=0}[x^{b_ip^k}](1+x^{p^k})^{a^i}\pmod p

\dbinom{n}{m}=\prod\limits_{i=0}[x^{b_i}](1+x)^{a^i}\pmod p

\dbinom{n}{m}=\prod\limits_{i=0}\dbinom{a_i}{b_i}\pmod p

即 : 把 n,m 按照 p 进制拆位,每一位对应求组合数乘积即可。

#include<cstdio>
#define MaxN 100500
#define ll long long
using namespace std;
int mod;
ll powM(ll a,int t=mod-2){
  ll ret=1;
  while(t){
    if (t&1)ret=ret*a%mod;
    a=a*a%mod;t>>=1;
  }return ret;
}
ll fac[MaxN],ifac[MaxN];
ll sC(int n,int m){
  if (n<m)return 0;
  return fac[n]*ifac[m]%mod*ifac[n-m]%mod;
}
ll C(ll n,ll m){
  if (!m)return 1;
  return sC(n%mod,m%mod)*C(n/mod,m/mod)%mod;
}
void Init()
{
  fac[0]=1;
  for (int i=1;i<mod;i++)
    fac[i]=fac[i-1]*i%mod;
  ifac[mod-1]=powM(fac[mod-1]);
  for (int i=mod-1;i;i--)
    ifac[i-1]=ifac[i]*i%mod;
}
void solve()
{
  ll n,m;
  scanf("%lld%lld%d",&n,&m,&mod);
  m+=n;Init();
  printf("%lld\n",C(m,n));
}
int main()
{
  int T;scanf("%d",&T);
  while(T--)solve();
  return 0;
}

例题 P4345 [SHOI2015]超能粒子炮·改

题意 : 多次给出 n,k ,求 \sum\limits_{i=0}^k\dbinom{n}{i}\pmod {2333}

设 p=2333 ,显然这是个素数。

根据卢卡斯定理有 \dbinom{n}{m}=\dbinom{\lfloor n/p\rfloor}{\lfloor m/p\rfloor}\dbinom{n\bmod p}{m\bmod p}

ANS=\sum\limits_{i=0}^k\dbinom{n}{i}

=\sum\limits_{i=0}^k\dbinom{\lfloor n/p\rfloor}{\lfloor i/p\rfloor}\dbinom{n\bmod p}{i\bmod p}

前半部分是块状变化的,我们分别枚举\lfloor i/p\rfloor 和 i\bmod p。注意最后多出来的一个散块。

=\sum\limits_{i=0}^{\lfloor k/p\rfloor-1}\dbinom{\lfloor n/p\rfloor}{i}\sum\limits_{j=0}^{p-1}\dbinom{n\bmod p}{j}+\dbinom{\lfloor n/p\rfloor}{\lfloor k/p\rfloor}\sum\limits_{j=0}^{k\bmod p}\dbinom{n\bmod p}{j}

注意到\sum\limits_{i=0}^{\lfloor k/p\rfloor-1}\dbinom{\lfloor n/p\rfloor}{i},\quad\sum\limits_{j=0}^{p-1}\dbinom{n\bmod p}{j},\quad\sum\limits_{j=0}^{k\bmod p}\dbinom{n\bmod p}{j}均是子问题。

设 f(n,k)=\sum\limits_{i=0}^k\dbinom{n}{i}

则有 f(n,k)=f(\lfloor n/p\rfloor,\lfloor n/k\rfloor-1)f(n\bmod p,p-1)+\dbinom{\lfloor n/p\rfloor}{\lfloor k/p\rfloor}f(n\bmod p,k\bmod p)

先预处理 p\times p 的 f,这用组合数递推不难做到 O(p^2).

然后就是求解 \dbinom{\lfloor n/p\rfloor}{\lfloor k/p\rfloor} ,使用Lucas定理即可做到 O(\log_p n)

观察参数 n 的变化 : 不断除以 p。

然后由于 f(n,k) 在 k>n 时等于 f(n,n) ,当 n\leq p 时就已经达到边界,只需要递归 O(\log_p n) 次。

复杂度O(T\log_p^2n+p^2).

提交记录

周边 : Kummer定理

注意到 $n!$ 中 $p$ 的次数是 $\sum\limits_{i=0}\lfloor n/p^i\rfloor$ ,证明比较显然。根据 $\dbinom{n+m}{m}=\dfrac{(n+m)!}{n!m!}

则因子 p 的个数是 \sum\limits_{i=0}\lfloor (n+m)/p^i\rfloor-\lfloor n/p^i\rfloor-\lfloor m/p^i\rfloor

当某一位上 \lfloor n/p^i\rfloor 相当于在 p 进制表示下去掉后 i 位。

那么只有 n+m 在这里进位了才能逃脱被去掉的命运而恰好多1.

例题 : CF582D Number of Binominal Coefficients

题意 : 给出参数 lim,p,a ,求有多少个组合数 \dbinom{n}{k} 是 p^a 的倍数,且 0\leq k \leq n \leq lim.

根据 Kummer 定理,能将问题转化成这样的形式:

求两个数 x,y 使得 x+y\leq lim ,而且会在 p 进制下产生多于 a 个进位。

这看起来就很数位DP的样子。

设f[i][j][0/1][0/1]表示考虑了前i位,还需要产生j个进位,是否卡上界,这一位是否进位。

注意这里是不需要记录前导 0 的,而且由于进位会影响上一位,还要钦定是否进位。

转移的时候分几类讨论 : (下一位进/不进位)(这一位进/不进位),分别计算可能的数对个数即可。

复杂度O(\log^2lim)

中段综合测试 : P2480 [SDOI2010]古代猪文

题意 : 给出n,g,求g^{\small\sum\limits_{k|n}\binom{n}{k}},对999911659取模。

首先,模数是个质数,根据费马小定理,指数要对 mod-1 取模。

则对 999911658=2\times 3\times 4679\times 35617 取模。

观察到没有重复的素因子,这比较和善。

求出分别模 2,3,4679,35617 的结果,然后使用中国剩余定理合并即可(甚至不需要Exgcd)。

求组合数可以使用Lucas定理。

复杂度是O(d(n)\log p+T(2+3+4679+35617))

扩展欧拉定理(降幂塔)

上集我们已经证明了 : [a\perp m]\ \Longrightarrow\ a^{\varphi(m)}=1\pmod m

现在我们来看扩展形式, \color{blue}\text{定理:}

a^n=\begin{cases}a^{n\bmod\varphi(m)}&(a\perp m)\\a^n&(a\not\perp m,\ n< \varphi(m))\\a^{(n\bmod \varphi(m))+\varphi(m)}&(a\not\perp m,\ n\geq\varphi(m))\end{cases}\large\pmod m

第一种情况已被证明,第二种情况是显然的,现在我们来证明第三种情况。

说人话就是 : 当 a\not\perp m 时,第一个 \varphi(m) 不是循环,从第二个才开始进入。

不妨先考虑 a=p^k 的情况。

此时若有 p\not\perp m ,则 m 可以表示成 s*p^c ,且 s\perp p。

那么根据普通欧拉定理有 p^{\varphi(s)}=1\pmod s 。同时,由于 \varphi 函数的积性,有 \varphi(s)|\varphi(m)。

于是有 p^{\varphi(m)}=1\pmod s

两边同时乘以 p^c ,能够得到 p^{\varphi(m)+c}=p^c\pmod {m}

注意到 \varphi(m)\geq \varphi(p^c)\geq p^{c-1}(p-1)\geq c

所以有 p^c=p^{c+\varphi(m)}=p^{c\bmod\varphi(m)+\varphi(m)}\pmod m ,我们证明了对于 p^c 满足结论。

结论同时也对于 p^{kc} 成立,可以归纳证明。

若有 2c\geq n\geq \varphi(m) \geq c ,能推知 n-c\leq \varphi(m)

则 p^n=p^{c}p^{n-c}=p^{c\bmod\varphi(m)+\varphi(m)+(n-c)\bmod\varphi(m)}

有可能会出现 c|n ,此时已经证毕。否则有 c\bmod\varphi(m)+(n-c)\bmod\varphi(m)=n\bmod\varphi(m)

然后对 n 的上界归纳即可。

现在我们对任意的 p^n 证明了 p^n=p^{(n\bmod\varphi(m))+\varphi(m)}

然后使用唯一分解定理即可把结论扩展至全体正整数。

例题 : P4139 上帝与集合的正确用法

题意 : 求 2^{2^{2^{2^{\dots}}}}\bmod m , \ p\leq 10^7

设G(p)=2^{2^{2^{2^{\dots}}}}\bmod m

根据扩展欧拉定理能得到 G(p)=2^{G(\varphi(p))+\varphi(p)}

边界就是 2^{2^{2^{2^{\dots}}}}\bmod 1=0

取多少次 \varphi(p) 会得到 1 呢?

注意到欧拉函数的其中一个定义式 : \varphi(m)=m\prod\limits_{p|m}\frac{p-1}{p}

当 m 无素因子 2 时,一个奇素因子一定会产生因子 2 ,因为 p-1 是偶数。

当 m 有素因子 2 时,值至少减半。

综上,经过 O(\log p) 次递归就能得到 1。

线性筛出欧拉函数复杂度是 O(p+T\log^2p) 的。

提交记录

附加题 :
- P4861 按钮
若 \gcd(k,m)>1 ,显然无解。

否则满足经典欧拉定理,则有 k^{\varphi(m)}=1\pmod m ,但是 \varphi(m) 并不是答案,可能存在更小的循环节。

可能的循环节一定是 \varphi(m) 的约数,否则能导出矛盾。

于是只需要判定 \varphi(m) 的约数即可,复杂度为 O(\sqrt{m}\log m)。

[提交记录] ()
- P3934 Nephren Ruq Insania
考虑乘方降幂塔,在 O(\log p) 层之后模数就会变为 1,所以我们只需要暴力取出 l 后的若干项即可。

上一题中,由于 2^{2^{2^{2^{\dots}}}} 恒大于 \varphi(m) ,所以总可以使用第三种情况。

现在由于可能产生第二种情况,我们计算时需要按照如下原则:

如果算的的数大于 m 则替换成 (a\bmod m)+m,否则不变。

提交记录
- CF906D Power Tower : 上一题的静态版,但数据较强。
- P3747 [六省联考2017]相逢是问候
仍然考虑结论 : 降幂塔不超过 O(\log p) 层。某个位置被赋值 O(\log p) 次之后就必然变为 c^{c^{c^{\dots}}}\bmod p。

考虑把不再变化的一整个区间冻结,可以证明只会产生 O(n\log p) 次修改叶子的操作。

每次单点修改需要计算一次降幂塔,朴素实现复杂度是O(\log p^2)的,无法通过。

模数只有 O(\log p) 个,可以分别于处理 c 的光速幂。

总复杂度为O(\sqrt{p}\log p+m\log ^2p)。

提交记录
- P5240 Derivation
P4777 【模板】扩展中国剩余定理（EXCRT）

仍然是若干个方程组\begin{cases}x=c_1\pmod {p_1}\\x=c_2\pmod {p_2}\\...\\x=c_m\pmod {p_m}\end{cases}

求x的最小解,但不保证p_{1...m}互质。

仍然考虑两两合并\begin{cases}x=c_1\pmod{p_1}\\x=c_2\pmod {p_2}\end{cases}

原先我们要直接构造某个模数在另一个模数下的逆,现在似乎不太行了。

第一个方程的通解为c_1+kp_1,现在就是求一个 k 使得c_1+kp_1=c_2\pmod{p_2}

移项有kp_1=c_2-c_1\pmod{p_2},熟练的同学已经知道怎么做了。

等价于kp_1+yp_2=c_2-c_1

根据斐蜀定理,方程有解的充要条件是\gcd(p_1,p_2)|(c_2-c_1).

如果可解,先扩欧求出kp_1+yp_2=\gcd(p_1,p_2),然后乘以\dfrac{c_2-c_1}{\gcd(p_1,p_2)}即可。

注意,新的模数是lcm(p_1,p_2)。此外还要注意防止溢出。

复杂度仍然是O(m\log p).

#include<cstdio>
#define ll long long
using namespace std;
inline ll mul(ll a,ll b,ll m){
  ll d=((long double)a/m*b+0.5);
  ll r=a*b-d*m;
  return r<0?r+m:r;
}
ll gcd(ll a,ll b)
{return !b ? a : gcd(b,a%b);}
void exgcd(ll a,ll b,ll &x,ll &y){
  if (b==0){x=1;y=0;return ;}
  exgcd(b,a%b,y,x);y-=(a/b)*x;
}
int n;
ll c1,p1,c2,p2;
int main()
{
  scanf("%d",&n);
  c1=0;p1=1;
  for (int i=1;i<=n;i++){
    scanf("%lld%lld",&p2,&c2);
    ll pd=gcd(p1,p2),sp=p2/pd*p1,k,y;
    c2=(c2-c1%p2+p2)%p2;
    //if (c2%pd) No Sol;
    exgcd(p1,p2,k,y);
    k=mul(k,(c2/pd),p2);
    c1=(c1+mul(p1,k,sp))%sp;
    p1=sp;
  }printf("%lld",c1);
  return 0;
}

例题 : P4774 [NOI2018]屠龙勇士

发现不会做这道题是写本文的动力之一。

首先容易发现,屠龙的规则和顺序严格确定,那么如果我们能成功,每一轮使用的剑是相同的。拿multiset维护一下便可。

设杀第i条龙的剑攻击力为k_i.

击杀每条龙的条件很像取模,我们能列出式子k_ix=a_i\pmod {p_i}

问题在于,可能有a_i>p_i,此时可能在比0大但为p_i倍数的血量停下来,所以x满足上述方程并不是充要条件。

考虑对每条龙记录砍到负数的最小次数的最大值,最后x在通解中取大一点就好了。

接下来的问题就是解方程组\begin{cases}k_1x=c_1\pmod {p_1}\\k_2x=c_2\pmod {p_2}\\...\\k_mx=c_m\pmod {p_m}\end{cases}

考虑将kx=c\pmod p变为形如x=c\pmod p的经典形式。

先特判0 : 当p|k且p|c时方程恒成立,不考虑。当p|k但p\!\!\not|\ c时无解。

由于k,p不一定互质,我们不能直接通过逆元来转化。

不过,我们仍然能尝试求解kx+py=c,无解则原方程组无解。

否则由斐蜀定理得\gcd(k,p)|c,设d=\gcd(k,p)。

则有\dfrac{k}{d}x+\dfrac{p}{d}y=\dfrac{c}{d}.

我们对\dfrac{p}{d}取模即得到\dfrac{k}{d}x=\dfrac{c}{d}\pmod {\dfrac{p}{d}}

由于\dfrac{k}{d}\perp\dfrac{p}{d},我们就可以求逆元变为经典形式了。

剩下的就是一个扩展中国剩余定理。

#include<algorithm>
#include<cstdio>
#include<set>
#define Itor set<ll>::iterator
#define ll long long
#define MaxN 100500
using namespace std;
const int mod=998244353;
ll mul(ll a,ll b,ll m)
{return (((a>>20)*b%m<<20)+(a-(a>>20<<20))*b)%m;}
ll gcd(ll a,ll b)
{return !b ? a : gcd(b,a%b);}
ll lcm(ll a,ll b)
{return a/gcd(a,b)*b;}
void exgcd(ll a,ll b,ll &x,ll &y){
  if (!b){x=1;y=0;return ;}
  exgcd(b,a%b,y,x);y-=(a/b)*x;
}
ll inv(ll a,ll m){
  ll x,y;exgcd(a,m,x,y);
  return (x+m)%m;
}
bool fl;
void merge(ll &c1,ll &p1,ll c2,ll p2)
{
  ll c=(c2-c1%p2+p2)%p2,d=gcd(p1,p2);
  if (c%d){puts("-1");fl=1;return ;}
  ll k,y;exgcd(p1,p2,k,y);
  ll p=lcm(p1,p2);
  c1=(c1+mul(mul(k,p1,p),c/d,p))%p;
  p1=p;
}
multiset<ll> s;
ll a[MaxN],k[MaxN],p[MaxN],t[MaxN]
  ,mx,tc,tp;
void calc(ll k,ll c,ll p)
{
  mx=max(mx,(c-1)/k+1);
  if (k%p==0){
    if (c%p==0)return ;
    else {puts("-1");fl=1;return ;}
  }
  ll d=gcd(k,p);
  if (c%d){puts("-1");fl=1;return ;}
  k/=d;p/=d;c/=d;
  c=mul(c,inv(k,p),p);
  merge(tc,tp,c,p);
}
int n,m;
void solve()
{
  scanf("%d%d",&n,&m);
  for (int i=1;i<=n;i++)scanf("%lld",&a[i]);
  for (int i=1;i<=n;i++)scanf("%lld",&p[i]);
  for (int i=1;i<=n;i++)scanf("%lld",&t[i]);
  s.clear();
  for (int i=1;i<=m;i++){
    ll x;scanf("%lld",&x);
    s.insert(x);
  }
  for (int i=1;i<=n;i++){
    Itor it=s.upper_bound(a[i]);
    if (it!=s.begin())it--;
    k[i]=*it;s.erase(it);
    s.insert(t[i]);
  }
  mx=0;tp=1;tc=0;fl=0;
  for (int i=1;i<=n&&!fl;i++)
    calc(k[i],a[i],p[i]);
  if (fl)return ;
  ll tk=mx/tp;
  while(tk*tp+tc<mx)tk++;
  printf("%lld\n",tk*tp+tc);
}
int main()
{
  int T;scanf("%d",&T);
  while(T--)solve();
  return 0;
}

P4720 【模板】扩展卢卡斯

跟普通卢卡斯定理关系不大 (

我们把模数分解成\prod p_i^{c_i}的形式,对于每个p^c单独求解,然后使用普通中国剩余定理合并即可。

考虑\dbinom{n}{m}=\dfrac{n!}{m!(n-m)!}\pmod{p^c}

但这毫无作用,往往有 p|(n!) ,这样就无法求逆了。

考虑将p的幂次提尽,设v(n)为n中含有p的幂次。

则有\dfrac{n!}{p^{v(n)}}≠0\pmod{p^c},这样就可以求逆了。设r(n)=\dfrac{n!}{p^{v(n)}}

我们计算\dfrac{r(n)}{r(m)r(n-m)}p^{v(n)-v(m)-v(n-m)}即可得到\dbinom{n}{m}。

上文介绍了v(n)=\sum\limits_{i=0}\lfloor n/p^i\rfloor ,现在问题变为求r(n).

对于n!=1*2*3*...*n,我们可以把其中为p的倍数的项提出来。

变为\Big[1*2*...*(p-1)\Big]*p*\Big[(p+1)(p+2)...(p+p-1)\Big]*2p*...

对于中间的每一块,\prod_{i=1}^{p-1}(kp+i)\pmod {p^c},可以变为[0,p^c)中的某一段连乘积,维护阶乘及其逆元即可。

这样的段会有O(n/p)个,但是由于\bmod\ p^c,我们可以批量计算。

具体见代码,注意散块的边界情况。这样计算一次的复杂度是O(p^c)的。

对于那些为p的倍数的位置,共\lfloor n/p\rfloor个,统一除以p之后,又变成了1*2*3*...

于是就变为求解r(\lfloor n/p\rfloor).

预处理复杂度为O(\sum p_i^{c_i}),回答一次的复杂度为O(\log^2 p).

( 如果采用欧拉定理+光速幂可以达到O(\sqrt{p}+\sum p_i^{c_i})预处理O(\log p)回答 )

#include<algorithm>
#include<cstdio>
#define MaxN 1005000
#define ll long long
using namespace std;
void exgcd(ll a,ll b,ll &x,ll &y){
  if (b==0){x=1;y=0;return ;}
  exgcd(b,a%b,y,x);y-=(a/b)*x;
}
ll inv(ll a,ll m){
  ll x,y;exgcd(a,m,x,y);
  return (x%m+m)%m;
}
ll powM(ll a,ll t,int mod)
{
  ll ret=1;
  while(t){
    if (t&1)ret=ret*a%mod;
    a=a*a%mod;t>>=1;
  }return ret;
}
ll v(ll n,int p){
  ll ret=0;
  while(n)ret+=(n/=p);
  return ret;
}
ll sav[MaxN];
ll r(ll n,int p,int mod){
  if (n==0)return 1;
  return powM(sav[mod-1],n/mod,mod)*sav[n%mod]%mod*r(n/p,p,mod)%mod;
}
ll C(ll n,ll m,int p,int mod)
{
  int c=v(n,p)-v(m,p)-v(n-m,p);
  ll ans=powM(p,c,mod);
  if (!ans)return 0;
  sav[0]=1;
  for (int i=1;i<mod;i++)
    if (i%p==0)sav[i]=sav[i-1];
    else sav[i]=sav[i-1]*i%mod;
  ans=ans*r(n,p,mod)%mod;
  ans=ans*inv(r(m,p,mod),mod)%mod;
  ans=ans*inv(r(n-m,p,mod),mod)%mod;
  return ans;
}
ll ret,mp=1;
void add(ll c,ll p)
{
  ll sp=mp*p;
  ret=(ret*inv(p,mp)%sp*p+mp*inv(mp,p)%sp*c)%sp;
  mp=sp;
}
ll n,m;int p;
int main()
{
  scanf("%lld%lld%lld",&n,&m,&p);
  int d=2;
  while(d*d<=p){
    if (p%d==0){
      int mod=1;
      while(p%d==0){
        p/=d;
        mod*=d;
      }add(C(n,m,d,mod),mod);
    }d++;
  }if (p>1)add(C(n,m,p,p),p);
  printf("%lld",ret);
  return 0;
}

例题① : P2183 [国家集训队]礼物

题意 : 有n个不同的小球,投入m个盒子里面,钦定第i个盒子要投w_i个球,问方案数。

对 mod 取模,保证 mod\leq 10^9 且分解之后的每个 p^c\leq 10^5

先判掉无解的情况。设S[k]=\sum\limits_{i=1}^kw_i

不难发现题目叫我们求的是\prod\limits_{i=1}^m\dbinom{n-S[i-1]}{w_i}

也就是在固定模数下求解m次组合数,套用EXLucas即可。

理论复杂度为O(\sqrt{p}+\sum p^c+m\log^2p),但是O\Big(m(\sqrt{p}+\sum p^c+\log^2p)\Big)也能过我就偷懒了。

评测记录

例题② : P3726 [AH2017/HNOI2017]抛硬币

题解Link

同余系基本定理2

P1495 【模板】中国剩余定理(CRT)

P3807 【模板】卢卡斯定理

扩展欧拉定理(降幂塔)

P4777 【模板】扩展中国剩余定理（EXCRT）

P4720 【模板】扩展卢卡斯