关于数论
关于数论
数论东西很多又很杂，所以想要总结一下，有一些算法的百度百科讲得很清楚，所以我就直接给了个链接在这（其实是懒23333），方便自己复习吧。
欧几里得算法
辗转相除法求gcd与lcm
使用辗转相除算出gcd后，lcm可以直接通过gcd算出，但是注意求lcm的过程可能爆int，建议使用long long
inline int gcd(int a,int b){while(b^=a^=b^=a%=b);return a;}
inline long long lcm(int a,int b){return (long long)a/gcd(a,b)*b;}
Tip : gcd(a,b) * lcm(a,b)=a * b
扩展欧几里得算法
int x,y;
inline int exgcd(int a,int b,int &x,int &y)
{
if(!b)
{
x=1,y=0;
return a;
}
int gcd=exgcd(b,a%b,x,y);
int t=x;x=y,y=t-a/b*y;
return gcd;
}
应用：
- 求不定方程：ax+by=c
- 解模线性方程：ax≡m(mod b) 转化为ax+by=m,与不定方程求解方法相同。
- 求解乘法逆元
逆元
这里多说几句逆元：
定义：存在ax≡1(mod p)，则称x是a关于模p的乘法逆元
定理：a关于p的乘法逆元的充要条件是gcd(a,p)=1
用途：因为两数相除不满足模的性质，所以用乘这个数的逆元代替：(a/b)%p=(a * x)%p
Tip : mod 1的逆元就是1
例题HDU1576
在题目中,我们将B在模9973下的逆元设为x，求(A/B)%9973就可以化为(A×x)%9973，这样的话就可以通过模的性质变成A%9973×x%9973，而在题目中A%9973已经给出，那么只需要通过拓展欧几里得求出逆元之后就可以得到答案了
代码：
#include<cstdio>
#define ll long long
using namespace std;
template<typename TP>inline void read(TP&x)
{
x=0;int f=1;char c=getchar();
while(c>'9'||c<'0'){if(c=='-')f=-1;c=getchar();}
while(c>='0'&&c<='9'){x=(x<<1)+(x<<3)+(c^48);c=getchar();}
x*=f;
}
template<typename TP>inline void print(TP x)
{
if(x<0)x=-x,putchar('-');
if(x>=10)print(x/10);
putchar(x%10+'0');
}
ll t,a,b,ans;
ll x,y;
inline ll exgcd(ll a,ll b,ll &x,ll &y)
{
if(!b)
{
x=1,y=0;
return a;
}
ll gcd=exgcd(b,a%b,x,y);
ll t=x;x=y,y=t-a/b*y;
return gcd;
}
int main()
{
read(t);
while(t--)
{
read(a),read(b);
exgcd(b,9973,x,y);
x=((x%9973)+9973)%9973;
print(a*x%9973),putchar('\n');
}
}
多数情况用扩展欧几里得就可以求得逆元了，当然还有用欧拉函数的方法，这里有一个时间复杂度为O(n)的神奇做法
只有一句递推式：
A[i]=-(p/i)*A[p%i];//p为模数
证明详见大佬的博客
博客
素数筛
最朴素的算法这里就不给了，直接上最优秀的
欧拉筛
利用了每个合数必有一个最小素因子，每个合数仅被它的最小素因子筛去正好一次。
prime数组中的素数是递增的,当 i 能被prime[j]整除 ，那么就不用继续重复筛除。
因为i%prime[j]=0，我们设i=prime[j]x，所以 i * prime[j+1]就等于 x * prime[j] * prime[j+1]
所以 i * prime[j+1] 必定被 prime[j]x 筛掉
代码中体现在：
if(i%prime[j]==0)break;
就可以做到几乎是线性的时间复杂度就可以筛出素数
int prime[2000000],size;
bitset<2000000>flag;
inline void make_prime(int x)
{
for(int i=2;i<=x;++i)
{
if(!flag[i]) prime[++size]=i;
for(int j=1;j<=size&&prime[j]*i<=x;++j)
{
flag[prime[j]*i]=1;
if(!(i%prime[j])) break;
}
}
}
费马平方和定理
内容如下：
除2以外，任何素数都能分成4k+1与4k+3的形式，任何形为4k+1的素数都能表示为两个平方数之和，而4k+3形式质数则不能
代码：
#include<cmath>
#include<bitset>
#include<vector>
#include<cstdio>
#include<cstring>
#include<iostream>
#define ri register int
#define rc register char
using namespace std;
template<typename TP>inline void read(TP&x)
{
x=0;int f=1;char c=getchar();
while(c>'9'||c<'0'){if(c=='-')f=-1;c=getchar();}
while(c>='0'&&c<='9'){x=(x<<1)+(x<<3)+(c^48);c=getchar();}
x*=f;
}
template<typename TP>inline void write(TP x)
{
if(x==0){putchar('0');return;}
char f[100]={0};int i=0;
if(x<0){putchar('-'),x=-x;}
while(x){++i,f[i]=x%10+48;x/=10;}
while(i){putchar(f[i]),--i;}
}
int prime[20000000],size;
bitset<300000005>flag;
int l,r,tot;
inline void make_prime(int x)
{
for(register int i=2;i<=x;++i)
{
if(!flag[i]) prime[++size]=i;
for(register int j=1;j<=size&&i*prime[j]<=x;++j)
{
flag[prime[j]*i]=1;
if(!(i%prime[j])) break;
}
}
}
int main()
{
read(l),read(r);
make_prime(r);
for(register int i=1;i<=size;++i)
{
if(prime[i]<l)continue;
if(prime[i]%4==1)
++tot;
}
if(l<=2&&r>=2) ++tot;
write(tot);
return 0;
}
欧拉函数
1到n中与n互质的数的个数被称为欧拉函数，记为φ(n)
1.通式：φ(x)=x(1-1/p1)(1-1/p2)(1-1/p3)(1-1/p4)…..(1-1/pn)
2.其中p1, p2……pn为x的所有质因数，x是不为0的整数。
3.特殊：φ(1)=1（唯一和1互质的数(小于等于1)就是1本身）
4.证明（容斥原理）：
p1,p2,p3...pk为n的质因子。
与n不互质的数的个数为:
n/p1+n/p2+...+n/pk-n/(p1 * p2)-...-n/(pk-1 * pk)-n/(p1 * p2 * p3)-...-n/(pk-2 * pk-1 * pk)-... +n/(p1 * p2 * ... * pk)
所以与n互质的数的个数为：
φ(n)=n-[n/p1+n/p2+...+n/pk-n/(p1 * p2)-...-n/(pk-1 * pk)-n/(p1 * p2 * p3)-...-n/(pk-2 * pk-1 * pk)-... +n/(p1 * p2 * ... * pk)]=n * (1-1/p1) * (1-1/p2) * ... * (1-1/pk);
性质
1.n>1时，1到n中与n互质的数的和为 n * φ(n)/2
2.若a，b互质，则φ(a * b)=φ(a) * φ(b)
3.若p为质数，则φ(p)=p-1(费马小定理)
4.若p为质数，p|n且p^2|n,则φ(n)=φ(n/p) * p
5.若p为质数，p|n且p^2不整除n,则φ(n)=φ(n/p) * (p-1)
求欧拉函数
求单个欧拉函数的方法
inline int eular(int n)
{
int ret=n;
for(int i=2;i<=sqrt(n);i++)
{
if(n%i==0)
{
ret=ret/i*(i-1);
while(n%i==0)n/=i;
}
}
if(n>1)ret=ret/n*(n-1);
return ret;
}
线性筛法
比较简单的证明可以通过上面的性质2、3、4、5推导，这里就不细说
int prime[2000000],eular[2000000],size;
bitset<2000000>flag;
inline void Eular(int x)
{
eular[1]=1;
for(int i=2;i<=x;++i)
{
if(!flag[i]) prime[++size]=i,eular[i]=i-1;
for(int j=1;j<=size&&prime[j]*i<=x;++j)
{
flag[prime[j]*i]=1;
if(!(i%prime[j]))
{
eular[prime[j]*i]=eular[i]*prime[j];
break;
}
eular[prime[j]*i]=eular[i]*eular[prime[j]];
}
}
}
用途
这里有一个重要的用途：用它来求一个数的逆元(首先要求这个数满足存在逆元)
在模为素数p的情况下，有费马小定理（欧拉定理在素数下的特殊情况）求得逆元
a^(p-1)≡1(mod p)
那么a^(p-2)*a≡1(mod p)
所以a的逆元为a^(p-2)
而在模不为素数的情况下，用欧拉定理求得逆元
a^φ(m)≡1(mod m)
同理 a^(φ(m)-1)*a≡1(mod m)
所以a的逆元为a^(φ(m)-1)
快速幂与慢速乘法
这两种算法都运用到了一些分治或者说是二进制的思想，并且都为了防止溢出而运用到了模的定理进行步步求模
快速幂
inline long long Pow(long long x,long long y)
{
long long ans=1;
for(long long i=y;i;i>>=1)
{
if(i&1)ans=ans*x%mod;
x=x*x%mod;
}
return ans%mod;
}
慢速乘法
inline long long Mul(long long x,long long y)
{
long long ans=0;
for(long long i=y;i;i>>=1)
{
if(i&1) ans=(ans+x)%mod;
x=(x+x)%mod;
}
return ans%mod;
}
还有一种比较暴力的溢出式取模，我不知道稳定性怎样，但是测试过后都是对的，并且速度会快很多233333
inline long long Mul(long long a,long long b,long long mod)
{
a%=mod,b%=mod;
long long c=(long double)a*b/mod;
long long ret=a*b-c*mod;
if(ret<0) ret+=mod;
else if(ret>=mod) ret-=mod;
return ret;
}
可以将快速幂与慢速乘法结合，解决更大的数据
inline long long Pow(long long x,long long y)
{
long long ans=1;
for(long long i=y;i;i>>=1)
{
if(i&1)ans=Mul(ans,x)%mod;
x=Mul(x,x)%mod;
}
return ans%mod;
}
组合数学
错排公式
错排公式：F[n]=(n-1)*(F[n-1]+F[n-2])
证明：
当n个编号元素放在n个编号位置，元素编号与位置编号各不对应的方法数用D(n)表示，那么D(n-1)就表示n-1个编号元素放在n-1个编号位置，各不对应的方法数，其它类推.
第一步，把第n个元素放在一个位置，比如位置k，一共有n-1种方法；
第二步，放编号为k的元素，这时有两种情况：⑴把它放到位置n，那么，对于剩下的n-1个元素，由于第k个元素放到了位置n，剩下n-2个元素就有D(n-2)种方法；⑵第k个元素不把它放到位置n，这时，对于这n-1个元素，有D(n-1)种方法；
综上得到
F(n) = (n-1)[F(n-2)+F(n-1)]
特殊地，F(1) = 0, F(2) = 1.
模板题：HDU1465
代码：
#include<cstdio>
#define ll long long
#define ri register int
using namespace std;
ll a,hhh[21];
int main()
{
hhh[1]=0,hhh[2]=1,hhh[3]=2;
for(ri i=4;i<=20;++i)
hhh[i]=(i-1)*(hhh[i-1]+hhh[i-2]);
while(~scanf("%lld",&a))
printf("%lld\n",hhh[a]);
return 0;
}
中国剩余定理
模板
LL prime[2333],remain[2333];
LL x,y;
inline LL exgcd(LL a,LL b,LL &x,LL &y)
{
if(!b)
{
x=1,y=0;
return a;
}
LL gcd=exgcd(b,a%b,x,y);
LL t=x;x=y,y=t-a/b*y;
return gcd;
}
inline LL crt()
{
LL m=1,ans=0;
for(ri i=1;i<=3;++i)
m*=prime[i];
for(ri i=1;i<=3;++i)
{
LL temp=m/prime[i];
exgcd(prime[i],temp,x,y);
ans=(ans+y*temp*remain[i])%m;
}
return (ans+m)%m;
}
扩展中国剩余定理
代码：
#include<cstdio>
#define ll long long
#define ri register int
#define rll register long long
template<typename TP>inline void read(TP&x)
{
x=0;int f=1;char c=getchar();
while(c<'0'||c>'9'){if(c=='-')f=-1;c=getchar();}
while(c>='0'&&c<='9'){x=(x<<3)+(x<<1)+(c^48);c=getchar();}
x*=f;
}
//********快读快输********
inline ll Mul(ll x,ll y,ll mod)//慢速乘法避免溢出
{
ll ans=0;
for(rll i=y;i;i>>=1)
{
if(i&1)ans=(ans+x)%mod;
x=(x+x)%mod;
}
return ans%mod;
}
ll n,x,y,ans,m;
ll prime[233333],remain[233333];
inline ll exgcd(ll a,ll b,ll &x,ll &y)//扩展欧拉定理
{
if(!b)
{
x=1,y=0;
return a;
}
ll gcd=exgcd(b,a%b,x,y);
ll t=x;x=y,y=t-a/b*y;
return gcd;
}
inline void excrt()
{
m=prime[1],ans=remain[1];
for(ri i=2;i<=n;++i)
{
ll temp=((remain[i]-ans)%prime[i]+prime[i])%prime[i];
ll gcd=exgcd(m,prime[i],x,y);
x=Mul(x,temp/gcd,prime[i]);
ans+=m*x;
m*=prime[i]/gcd;
ans=(ans+m)%m;
}
printf("%lld",ans);
}
int main()
{
read(n);
for(ri i=1;i<=n;++i)
read(prime[i]),read(remain[i]);//输入模数与余数
excrt();
return 0;
}
这里是关于数论的洛谷日报整理：
