数论是我一生之敌。QAQ
这个巨佬的博客数论是极好的,但是,他好像弃坑了。
欧拉筛
emmm,这个比较简单,就贴代码。
for (i=2; i<maxn; ++i) {
if (!isn_prime[i]) prime[++cnt]=i;
for (phi[i]=cnt,j=1; j<=cnt&&prime[j]*i<maxn; ++j) {
isn_prime[i*prime[j]]=1;
if (i%prime[j]==0) break;
}
}
扩欧求逆元
推导
老是会忘掉。qwq
a x + b y = g c d ( a , b ) = 1 ax+by=gcd(a,b)=1 ax+by=gcd(a,b)=1
∵ g c d ( a , b ) = g c d ( b , a % b ) \because gcd(a,b)=gcd(b,a\%b) ∵gcd(a,b)=gcd(b,a%b)
∴ a x + b y = g c d ( a , b ) = g c d ( b , a % b ) = b x ′ + ( a % b ) y ′ \therefore ax+by=gcd(a,b)=gcd(b,a\%b)=bx'+(a\%b)y' ∴ax+by=gcd(a,b)=gcd(b,a%b)=bx′+(a%b)y′
⇒ b x ′ + ( a − b ∗ ⌊ a b ⌋ ) y ′ \Rightarrow bx'+(a-b*\lfloor\cfrac{a}{b}\rfloor)y' ⇒bx′+(a−b∗⌊ba⌋)y′
⇒ a y ′ + b ( x ′ − y ′ ∗ ⌊ a b ⌋ ) \Rightarrow ay'+b(x'-y'*\lfloor\cfrac{a}{b}\rfloor) ⇒ay′+b(x′−y′∗⌊ba⌋)
∴ x = y ′ , y = ( x ′ − y ′ ∗ ⌊ a b ⌋ ) \therefore x=y',y=(x'-y'*\lfloor\cfrac{a}{b}\rfloor) ∴x=y′,y=(x′−y′∗⌊ba⌋)
代码
int ex_gcd(int a,int b,int &x,int &y) {
if (b==0) return x=1,y=0,a;
int ans=ex_gcd(b,a%b,x,y);
int tmp=x;
x=y,y=tmp-a/b*y;
return ans;
}
线性求逆元
i n v [ a ] = i n v [ P % a ] ∗ ( P − ⌊ P a ⌋ ) inv[a]=inv[P\% a]*(P-\lfloor\cfrac{P}{a}\rfloor) inv[a]=inv[P%a]∗(P−⌊aP⌋)
证明
令 x = P % a , y = ⌊ P a ⌋ x=P\% a,y=\lfloor\cfrac{P}{a}\rfloor x=P%a,y=⌊aP⌋
则 x + a y = P x+ay=P x+ay=P
即 x + a y ≡ 0 ( m o d P ) x+ay\equiv0\ \ (mod\ P) x+ay≡0 (mod P)
x ≡ − a y ( m o d P ) x\equiv-ay\ (mod\ P) x≡−ay (mod P)
i n v ( a ) ≡ i n v ( x ) ∗ ( P − y ) ( m o d P ) inv(a)\equiv inv(x)*(P-y)\ (mod\ P) inv(a)≡inv(x)∗(P−y) (mod P)
∴ i n v ( a ) ≡ i n v ( P % a ) ∗ ( P − ⌊ P a ⌋ ) \therefore inv(a)\equiv inv(P\%a)*(P-\lfloor\cfrac{P}{a}\rfloor) ∴inv(a)≡inv(P%a)∗(P−⌊aP⌋)
代码
for (inv[1]=1,i=2; i<maxn&&i<P; ++i) inv[i]=1ll*inv[P%i]*(P-P/i)%P;
例题
SDOI2008 沙拉公主的困惑
思路
我好弱,FTQ说这种题很经典,我却没做过。QwQ
我们可以吧 N ! N! N!分成长 M ! M! M!块,显然每块答案是 ϕ ( M ! ) \phi(M!) ϕ(M!)
所以 a n s = N ! M ! ∗ φ ( M ! ) ans=\cfrac{N!}{M!}*\varphi(M!) ans=M!N!∗φ(M!)
可以知道对于 φ ( N ) \varphi(N) φ(N)有
φ ( N ) = φ ( ∏ 1 m P i ε i ) = ∏ i m ( P i ε i ∗ P i − 1 P i ) = N ∗ ∏ i m ( P i − 1 P i ) \varphi(N)=\varphi(\prod_{1}^{m}P_i^{\varepsilon_{i}})=\prod_i^m(P_i^{\varepsilon_i}*\cfrac{P_i-1}{P_i})=N*\prod_i^m(\cfrac{P_i-1}{P_i}) φ(N)=φ(∏1mPiεi)=∏im(Piεi∗PiPi−1)=N∗∏im(PiPi−1)
所以 a n s = N ! ∗ ∏ ( P − 1 ) ∗ ∏ i n v ( P ) ans=N!*\prod(P-1)*\prod inv(P) ans=N!∗∏(P−1)∗∏inv(P) P为小于等于 M M M的质数。
代码
#include<cstdio>
#include<string>
const int maxn=1e7+5;
bool isn_prime[maxn];
int N,M,T,P,cnt;
int prime[maxn],phi[maxn];
int inv[maxn],fac[maxn];
int fac_prime[maxn],fac_inv_prime[maxn];
void init() {
register int i,j;
for (i=2; i<maxn; ++i) {
if (!isn_prime[i]) prime[++cnt]=i;
for (phi[i]=cnt,j=1; j<=cnt&&prime[j]*i<maxn; ++j) {
isn_prime[i*prime[j]]=1;
if (i%prime[j]==0) break;
}
}
for (inv[1]=1,i=2; i<maxn&&i<P; ++i) inv[i]=1ll*inv[P%i]*(P-P/i)%P;
for (fac_prime[0]=i=1; i<=cnt; ++i)
fac_prime[i]=1ll*fac_prime[i-1]*(prime[i]-1)%P;
for (fac_inv_prime[0]=i=1; i<=cnt; ++i)
if (prime[i]==P) fac_inv_prime[i]=fac_inv_prime[i-1];
else fac_inv_prime[i]=1ll*fac_inv_prime[i-1]*inv[prime[i]%P]%P;
for (fac[0]=i=1; i<maxn; ++i)
if (i!=P) fac[i]=1ll*fac[i-1]*i%P;else fac[i]=fac[i-1];
}
int main() {
scanf("%d%d",&T,&P),init();
while (T--) {
scanf("%d%d",&N,&M);
if (N>=P&&M<P) puts("0");else
printf("%lld\n",1ll*fac[N]*fac_prime[phi[M]]%P*fac_inv_prime[phi[M]]%P);
}
return 0;
}