QBXT游记 | Day2 Moring
Miller Rabin判断质数
首先需要知道的一点是这个方法的正确性并不是一定的，也就是说在一定的数据下，这个算法是会出错的
不过只要合理的选取质数就能在longlong范围内做到正确，具体的在代码里面加了注释
先说一下下面要用到的一个性质，也就是此法判断质数的核心
若n是质数，取\(1≤a<n\)，则如下的两条性质有至少一条成立：
$1. a^d \equiv 1 \pmod n $
\(2.a^{d×2^i}\equiv n-1\pmod n\)
钟长者说这是扩展费马小定理
long long ksc(long long x,long long y,long long n)
//计算x*y%n，实际上这并不是快速乘，只不过是能够避免爆int，这实际上是比朴素的乘法要慢的
{
if (y==0) return 0;
long long z=ksc(x,y/2,n);//
z=(z+z)%n;
if (y%2==1) z=(z+x)%n;
return z;
}
long long ksm(long long x,long long y,long long n)
{
if(!y) Heriko 1%n;
long long e=ksm(x,y/2,n);
long long ans=e*e%n;
if(ans%2==1) ans=ans*x%n;
Heriko ans;
}
bool miller_rabin(long long n,int a)
//检查a这个数 能否满足两个性质中的至少一个
{
long long d=n-1;
int r=0;
while (d%2==0)
{
d=d/2;
r++;
}
long long x=ksm(a,d,n);//a^d % n
if (x==1) return true;
for (int i=0;i<r;i++)
{
if (x == n-1) return true;
x = ksc(x,x,n);
}
return false;
}
//O(k*logn*logn)
int prime_list[] = {2,3,5,11,37,67,73,97};
bool prime(long long n)
{
if (n<2) return false;
for (int a=0;a<8;a++)
{
if (n == prime_list[a]) return true;
if (n%prime_list[a] ==0) return false;
if (miller_rabin(n,prime_list[a]) == false) return false;
}
return true;
}
当然还有一些奇奇怪怪的筛法，一并放在这里
bool prime(int x)
{
if (x<2) return false;
for (int a=2;a*a<=x;a++)
if (x%a==0) return false;
return true;
}
for (int a=2;a<=n;a++)
for (int b=a+a;b<=n;b+=a)
not_prime[b] = true;//O(nlogn)
for (int a=2;a<=n;a++)
if (!not_prime[a])
for (int b=a+a;b<=n;b+=a)
not_prime[b] = true;//O(nloglogn)
Exgcd
扩展欧几里得算法 Or 辗转相除法 至尊纪念版
一般来说和解决二元一次不定方程和同余方程脱不了干系
下面贴上洛谷的此例题的代码（和题解可能相像)
#include <bits/stdc++.h>
#define Heriko return
#define Deltana 0
#define S signed
#define U unsigned
#define LL long long
#define R register
#define I inline
#define D double
#define LD long double
#define mst(a,b) memset(a,b,sizeof(a))
#define ON std::ios::sync_with_stdio(false)
using namespace std;
I void fr(LL &x)
{
LL f=1;char c=getchar();
x=0;
while(c<'0'||c>'9')
{
if(c=='-') f=-1;
c=getchar();
}
while(c>='0'&&c<='9')
{
x=(x<<3)+(x<<1)+c-'0';
c=getchar();
}
x*=f;
}
void fp(LL x)
{
if(x<0) x=-x,putchar('-');
if(x>9) fp(x/10);
putchar('0'+x%10);
}
LL gcd(LL x,LL y)
{
Heriko !y ? x : gcd(y,x%y);
}
void Exgcd(LL &x,LL &y,LL a,LL b)
{
if(!b)
{
x=1;y=0;
Heriko;
}
LL p;
Exgcd(x,y,b,a%b);
p=x;x=y;y=p-(a/b)*y;
Heriko;
}
LL t;
S main()
{
fr(t);
while(t--)
{
LL x,y,minx,miny,a,b,c,g,maxx,maxy,sum;
x=y=a=b=c=g=maxx=maxy=sum=0;minx=miny=0x3f;
fr(a),fr(b),fr(c);g=gcd(a,b);
if(c%g!=0) printf("-1\n");
else
{
a/=g,b/=g,c/=g;Exgcd(x,y,a,b);x*=c,y*=c;
minx= x>0 && x%b!=0 ? x%b : x%b+b;
maxy=(c-minx*a)/b;
miny= y>0 && y%a!=0 ? y%a : y%a+a;
maxx=(c-miny*b)/a;
if(maxx>0) sum=(maxx-minx)/b+1;
if(!sum) printf("%lld %lld\n",minx,miny);
else printf("%lld %lld %lld %lld %lld\n",sum,minx,miny,maxx,maxy);
}
}
Heriko Deltana;
}
以及上课时钟长者的Code
int ex_gcd(int a,int b,int &x,int &y)
{
if (b==0)
{
x=1;y=0;
return a;
}
int xp,yp;
int g=ex_gcd(b,a%b,xp,yp);
x = yp;
y = xp - yp * (a/b);
return g;
}
实际上这和我的原来的写法不一样，下面放上我之前的写法
I void Exgcd(LL a,LL b,LL &d,LL &x,LL &y)
{
if(b==0) d=a,x=1,y=0;
else Exgcd(b,a%b,d,y,x),y-=x*(a/b);
}
当然以上的实现方式大同小异，可以放心食用
在求同余方程的特解的时候有一种钟长者独创的据说除了他其他人都不会卡的大数翻倍法
下面贴出Code
void Merrge(int &a1,int &m1,int a2,int m2)
{
if(m2>m1) swap(m1,m2),swap(a1,a2);//大数翻倍法
while(a1%m2!=a2) a1+=m1;
m1=lcm(m1,m2);
}
这里之所以函数名字时Merrge是因为我编译的时候用的万能头，而algorithm里面merge是并归排序的STL函数
BSGS算法
一种求解\(a^x\equiv b\pmod p\)的算法，详情可以查阅OI WIKI中有关BSGS的介绍
Oi wiki 上面给的翻译是大步小步，钟长者说这是北上广深（
大概思路就是判断要找的数在那一行
//以下By Deltana
#include <bits/stdc++.h>
#define Heriko return
#define Deltana 0
#define S signed
#define U unsigned
#define LL long long
#define R register
#define I inline
#define D double
#define LD long double
#define mst(a,b) memset(a,b,sizeof(a))
#define ON std::ios::sync_with_stdio(false)
using namespace std;
I void fr(LL &x)
{
LL f=1;char c=getchar();
x=0;
while(c<'0'||c>'9')
{
if(c=='-') f=-1;
c=getchar();
}
while(c>='0'&&c<='9')
{
x=(x<<3)+(x<<1)+c-'0';
c=getchar();
}
x*=f;
}
int ksm(int x,int y,int n)
{
if(!y) Heriko 1%n;
int e=ksm(x,y/2,n);
int ans=e*e%n;
if(ans%2==1) ans=ans*x%n;
Heriko ans;
}
int bsgs(int a,int b,int p)
//求解a^x %p =b
{
int s=sqrt(p);
set<int> se;
int x=1;
for(R int i=0;i<s;i++)
{
se.insert(x);
x=1ll*x*a%p;
}
int y=ksm(ksm(a,s,p),p-2,p);//(a^s)^(-1)
int z=b;
for(R int i=1;;i++)
{
if(se.count(z))
{
z=ksm(a,(i-1)*s,p);
for(R int j=(i-1)*s;;j++)
{
if(z==b) Heriko j;
z=1ll*z*a%p;
}
}
z=1ll*z*y%p;
}
}
S main()
{
Heriko Deltana;
}
//以下By 钟长者
#include<cmath>
#include<set>
int bsgs(int a,int b,int p)
//求解 a^x % p =b
{
int s=sqrt(p);
set<int> se;
int x=1;
for (int i=0;i<s;i++)
{
se.insert(x);
x=1ll*x*a%p;
}
int y = ksm(ksm(a,s,p) , p-2,p);//(a^s)^(-1)
int z=b;
for (int i=1;;i++)
{
if (se.count(z) != 0)
{
z = ksm(a,(i-1)*s,p);
for (int j=(i-1)*s;;j++)
{
if (z == b) return j;
z = 1ll * z *a %p;
}
}
z = 1ll *z *y %p;
}
}
Do you like WHAT YOU SEE ?
