[学习笔记]数论基础篇2

原根

已知当$gcd(a,b)=1$时，$a^{\varphi (b)}\%b=1$
定义在模$b$意义下，$a$的阶为$d$，且$a^d\%b=1$
如果$a$的阶为$\varphi (b)$，则$a$就为$b$的原根

大步小步算法(bsgs)

$bsgs$算法用于求离散对数，即已知$a,b,p$求满足$a^x\equiv b(modp)$的$x$,其中$gcd(a,p)=1$
由欧拉定理可知$a^{\varphi (p)}\%p=1$,所以我们可以通过暴力枚举$x$来求是否满足条件，时间复杂度为$O(p)$
当$p$很大时，这个方法就直接被干掉了，于是有了$bsgs$算法
令$m=\lceil \sqrt{}p\rceil$,设$r=x\%m$ 则$x=k\ast m+r$ 其中$0<=k<m0<=r<m$
我们便有了$a^{k\ast m+r}\equiv b(modp)$
得到$a^{k\ast m}\equiv b\ast a^{-r}(modp)$
那么我们就可以处理出$b\ast a^{-r}\%p$的值，存进$map$里面
在枚举$k$，求到$a^{k\ast m}\%p$是否存在，存在则就可以输出$x=k\ast m+r$，如果都没有，则无解
在$a^{k\ast m}\equiv b\ast a^{-r}$中我们涉及到了要处理逆元，但是必须是要$gcd(a^r,p)=1$，但是一些丧心病狂的出题人就会卡你,所以我们设$x=k\ast m-r$ 其中$0<k<=m0<r<=m$
上面的式子就可以变成$a^{k\ast m}\equiv b\ast a^r(modp)$
处理便一样，记得要特判掉$b=0$的情况
代码是网上的一个巨佬的，因为找不到自己的了

LL m=ceil(sqrt(c));// 注意要向上取整
        mp.clear();
        if(a%c==0)
        {
        printf("no solution\n");
        continue;
        }
        // 费马小定理的有解条件
        LL ans;//储存每一次枚举的结果 b* a^j
        for(LL j=0;j<=m;j++)  // a^(i*m) = b * a^j
        {
            if(j==0)
            {
                ans=b%c;
                mp[ans]=j;// 处理 a^0 = 1
                continue;
            }
            ans=(ans*a)%c;// a^j
            mp[ans]=j;// 储存每一次枚举的结果
        }
        LL t=fastpow(a,m,c);
        ans=1;//a ^(i*m)
        LL flag=0;
        for(LL i=1;i<=m;i++)
        {
            ans=(ans*t)%c;
            if(mp[ans])
            {
                LL out=i*m-mp[ans];// x= i*m-j
                printf("%lld\n",(out%c+c)%c);
                flag=1;
                break;
            }

        }
        if(!flag)
        printf("no solution\n");

扩展大步小步算法（exbsgs）

扩展大步小步也是解决 已知$a,b,p$求满足$a^x\equiv b(modp)$的$x$这一类问题的，只是当$a,p$不互质的情况，上面的方法解决不了了。
根据给出的式子我们可以得到$a^x+kp=b$
因为$a,p$不互质，我们设$d=gcd(a,p)$
上式同除$d$,得到$a^{x-1}\ast a/d+kp/d=b/d$
当$b/d\ne 0$时，这个方程便无解
若$a,p^{\prime }$仍然不互质，继续除去$gcd(a,p^{\prime })$
当他们互质的时候，我们便有了$a^{x-cnt}\equiv b/a^{\prime }(mod{p}^{\prime })$
这里的$cnt$表示总共取了多少个$a$出来取$gcd$,$a^{\prime }$表示这些取出来的$a$除了$gcd$后的乘积
这里$a,p^{\prime }$又互质了，便回到了$bsgs$，为了避免逆元，我们便把左边的初始值设为$a^{\prime }$。

int qkpow(int a,int b){
    int ans=1;
    while(b){
        if(b&1)
            ans=1ll*ans*a%p;
        b>>=1;
        a=1ll*a*a%p;
    }
    return ans;
}
map<LL,int> mp;
void exbsgs(int x,int y){
    if(y==1){
        printf("0\n");
        return;
    }
    int d=__gcd(x,p),k=1,cnt=0;
    while(d!=1){
        if(y%d!=0){
            puts("No Solution");
            return;
        }
        cnt++;y/=d,p/=d;
        k=1ll*k*x/d%p;
        if(y==k){
            writeint(cnt);
            puts("");
            return;
        }
        d=__gcd(x,p);
    }
    int s=y;
    mp.clear();
    int m=sqrt(p-1)+1;
    for(int i=0;i<=m;i++){
        mp[s]=i;
        s=1ll*s*x%p;
    }
    s=k,k=qkpow(x,m);
    for(int i=1;i<=m;i++){
        s=1ll*s*k%p;
        if(mp.count(s)){
            writeint(i*m-mp[s]+cnt);
            puts("");
            return;
        }
    }
    puts("No Solution");
}

miller_rabbin算法

这是一个用来判断素数的方法，利用了随机算法可以在很快的时间内判断一个数是否是质数。
根据费马小定理，如果$p$是一个质数，则有$a^{p-1}\%p=1a<p$
但是反过来就不一样了，因为有一类称为迈克尔数也满足这个式子，但这个数不一定是质数
我们令$p-1=2^k\ast r$ $a^{p-1}=a^{r^{2k}}$
我们初始值设为$a^r$,然后不断平方
如果第一次出现模$p$的值为$1$，那么检测上一次的值是否为$-1$ ,如果是,则没有毛病；否则可以断定不是质数。
设$y^2\%p=1$
则$y^2-1\%p=0$
则$(y-1)\ast (y+1)\%p=0$
因为$gcd(a,p)=1$,所以$gcd(y,p)=1$
所以要么$(y-1)\%p=0$
要么$(y+1)\%p=0$
$y=1$或者$y=p-1$(即在在模$p$的意义下为$-1$)

#include <cstdio>
#include <cstdlib>
#include <ctime>
using namespace std;
#define int long long
int read()
{
    int num=0,flag=1;char c;
    while((c=getchar())<'0'||c>'9')if(c=='-')flag=-1;
    while(c>='0'&&c<='9')num=(num<<3)+(num<<1)+(c^48),c=getchar();
    return num*flag;
}
int n;
int ksc(int a,int b)
{
	int res=0;
	while(b)
	{
		if(b&1) res=res+a;
		a=a+a;
		if(res>=n) res-=n;
		if(a>=n) a-=n;
		b>>=1;
	}
	return res;
}
int qkpow(int a,int b)
{
	int r=1;
	while(b>0)
	{
		if(b&1) r=ksc(r,a);
		a=ksc(a,a);
		b>>=1;
	}
	return r;
}
int miller()
{
	int k=n-1,t=0;
	while(k%2==0)
	{
		t++;
		k/=2;
	}
	for(int i=0;i<=10;i++)
	{
		int a=rand()%(n-1)+1,ls,np=0;
		a=qkpow(a,k);ls=a;
		if(a==1) continue;
		for(int j=1;j<=t;j++)
		{
			a=ksc(a,a);
			if(a==1)
			{
				if(ls!=n-1) np=1;
				break;
			}
			ls=a;
		}
		if(np || a!=1) return 0;
	}
	return 1;
}
signed main()
{
	srand(time(0));
	while(~scanf("%I64d",&n))
	{
		if(n&1==0)
		{
			puts("It is not a prime number.");
			continue;
		}
		if(miller()) puts("It is a prime number.");
		else puts("It is not a prime number.");
	}
}

高斯消元

int gauss(){
	int row=1,col=1;
	for(;row<=n&&col<=m;row++,col++){
		int cur=row;
		for(int i=row+1;i<=n;i++){
			if(Fabs(a[i][col])>Fabs(a[cur][col])){
				cur=i;
			}
		}
		if(cur!=row){
			for(int i=col;i<=m+1;i++) 
				swap(a[cur][i],a[row][i]);
		}
		if(a[row][col]==0){
			row--;
			continue;
		}
		for(int i=row+1;i<=n;i++){
			if(a[i][col]!=0){
				LL LCM=lcm(a[i][col],a[row][col]);
				LL tmpa=LCM/a[i][col];
				LL tmpb=LCM/a[row][col];
				for(int j=col;j<=m+1;j++) 
					a[i][j]=(a[i][j]*tmpa-a[row][j]*tmpb)%mo;
			}
		}
	}
	for(int i=row;i<=n;i++) 
		if(a[i][col]!=0) 
			return -1;
	if(row<=m) return m-row+1;
	for(int i=n;i>=1;i--){
		LL tmp=a[i][m+1];
		for(int j=i+1;j<=m;j++) 
			tmp=(tmp-x[j]*a[i][j])%mo;
		x[i]=tmp*inv(a[i][i])%mo;
	}
	return 0;
}

#include <bits/stdc++.h>
using namespace std;
const int maxn = 105;
const double eps = 1e-7;
int readint(){
	int x=0,f=1;char s=getchar();
	#define sc (s=getchar())
	while(s<'0'||s>'9'){
		if(s=='-')
			f=-1;
		sc;
	}
	while(s>='0'&&s<='9'){
		x=(x<<3)+(x<<1)+(s^48);
		sc;
	}
	#undef sc
	return x*f;
} 
double a[maxn][maxn];
double ans[maxn];
int main (){
	int n=readint();
	for(int i=1;i<=n;i++){
		for(int j=1;j<=n+1;j++)
			scanf("%lf",&a[i][j]);
	} 
	for(int i=1;i<=n;i++){
		int now=i;
		for(int j=i+1;j<=n;j++){
			if(fabs(a[j][i])>fabs(a[now][i]))
				now=j;
		}
		if(fabs(a[now][i])<eps){
			puts("No Solution");
			return 0;
		}
		if(i!=now)
			swap(a[i],a[now]);
		for(int j=i+1;j<=n;j++){
			double temp=a[j][i]/a[i][i];
			for(int k=i;k<=n+1;k++)
				a[j][k]-=a[i][k]*temp;
		}
	}
	ans[n]=a[n][n+1]/a[n][n];
	for(int i=n-1;i>=1;i--){
		for(int j=1;j<=i;j++){
			a[j][n+1]-=ans[i+1]*a[j][i+1];
		}
		ans[i]=a[i][n+1]/a[i][i];
	}
	for(int i=1;i<=n;i++){
		printf("%.2f\n",ans[i]);
	}
	return 0;
}