zoukankan      html  css  js  c++  java
  • 基础数论板子整理

    前排预警

    本文几乎没有理论(也就是不能当做课件来看), 只有板子, 纯复习用。

    由于我很菜, 并且时间不太够, 于是只学了一部分。

    像什么 BSGS、ExBSGS、 Miller-Rabin, 还有那个大数质因数分解我都不会(淦)。

    一些基本操作都是现写的伪代码, 都标了(就是可能出错); 其他的Luogu有模板题的都是过了模板题的代码, 都标了

    有可能会写错

    前言

    基础数论的大部分内容是取模意义下的运算, 猜测最初引入这些科技的目的是为了代替高精度算法qwq。



    基本操作

    • (O(sqrt n))时间判断(n)是否是质数。

    理论支撑 : 若 (n) 为合数, 则其至少有一个小于等于(sqrt n)的质因子。

    代码():

    bool is_pm(int n) {
    	int tmp=sqrt(n);
    	for(int i=2; i<=tmp; ++i)
    		if(n%i == 0) return false;
    	return true;
    }
    
    • (O(sqrt n))时间将(n)分解质因数。

    理论支撑:(n)至多有一个大于(sqrt n)的质因子(不是一种)。

    代码():

    int tot=0, p[_], c[_]; //p[i]保存n的第i个质因子, c[i]保存p[i]的次数
    void divide(int n) {
    	tot = 0;
    	for(int i=2; i<=sqrt(n); ++i) if(n%i==0) {
    		//不用担心i是否为质数, 因为如果i是合数且i整除n
    		//那么i的质因子一定整除n, 所以枚举到合数时不会出现 n%i == 0
    		p[++tot] = i, c[tot] = 0;
    		while(n%i==0) ++c[tot], n/=i;
    	}
    	if(n>1) p[++tot] = n, c[tot] = 1;
    	return;
    }
    

    由于上面两个操作都是依赖于枚举质数的, 所以可以预处理一张质数表来加速枚举。

    由于([1,n])中的质数个数约为(frac{n}{ln(n)}), 所以加速后上面两个操作的时间复杂度都为(O(frac{sqrt n}{ln(sqrt n)}))(约为qwq)。


    • (O(log_2 a)) or (O(log_2 b)) 时间求(a、b)的最大公约数。(反正就是(O(log))级别的啦qwq)

    理论支撑就是代码里的式子(不会证)

    代码(

    int gcd(int a, int b) {
    	return (!b) ? a : gcd(b, a%b);
    }
    

    • (O(n))时间筛出([1,n])的素数表(还可以筛欧拉函数, 好像还可以筛莫比乌斯函数)

    理论支撑:将筛一个合数(x)的方式限定为用这个数的最小质因子(p)(frac{x}p)筛这个数, 意即将每个合数被筛的次数限定为一次, 故时间复杂度为(O(n))

    筛素数代码():

    int m, prime[maxn + 15], v[maxn + 15];
    //prime[i]表示第k大素数, v[i]表示数i的最小质因子
    void euler(int n) { //筛出[1,n]所有素数
    	for(int i=2; i<=n; ++i) {
    		if(!v[i]) v[prime[++m] = i] = i; //i为质数
    		for(int j=1; j<=m; ++j) {
    			if(prime[j] > n/i || prime[j] > v[i]) break;
    			//不筛大于n的, 不以限定形式之外的形式筛
    			v[prime[j] * i] = prime[j];
    		}
    	}
    	return;
    }
    

    exgcd

    • (O(log_2 a)) or (O(log_2 b)) 求解 (ax + by = gcd(a,b))的一组解(x、y)
      就是exgcd啦。

    理论支撑:是可以考场上推出来的东西所以不写了

    还有就是一些细节什么的做做题就可以回想起来了。
    青蛙的约会
    荒岛野人

    代码(

    void exgcd(int a, int b, int &x, int &y) {
    	//exgcd可以顺便求gcd(a,b), 且加上这个功能后代码还是一行, 但是我懒得写。
    	!b ? x=1,y=0 : (exgcd(b, a%b, y, x), y -= x*(a/b));
    }
    

    逆元相关

    (a)在膜(p)意义下有逆元的充要条件是(a perp p)

    • (O(log_2 a)) or (O(log_2 p))(a)在膜(p)意义下的逆元,条件是(aperp p)

    就是用exgcd解(ax equiv 1 (mod p))


    欧拉定理

    [a^{varphi(m)} equiv 1 ag{$a perp m$} ]


    • (O(sqrt n))时间求 (varphi(n))

    理论支撑:欧拉函数的公式。

    代码():

    int phi(int n) {
    	int res = n;
    	for(int i=2; i<=sqrt(n); ++i) if(n%i==0) {
    		res = res / i * (i-1);
    		while(n%i==0) n/=i;
    	}
    	if(n>1) res = res / i * (i-1);
    	return res;
    }
    

    • (O(n))([1,n])的欧拉函数
      理论支撑:欧拉函数的公式(自己推推)

    代码():

    int m, prime[maxn], v[maxn], phi[maxn];
    //phi[i]为i的欧拉函数
    //其他参考上面的筛素数代码
    void euler(int n) {
    	for(int i=2; i<=n; ++i) {
    		if(!v[i]) {
    			v[prime[++m] = i] = i;
    			phi[i] = i-1;
    		}
    		for(int j=1; j<=m; ++j) {
    			if(prime[j] > n/i || prime[j] > v[i]) break;
    			v[prime[j] * i] = prime[j];
    			phi[prime[j] * i] = phi[i] * (i%prime[j] ? prime[j]-1 : prime[j]);
    		}
    	}
    }
    

    • (O(n))([1,n-1])在膜(n)意义下的逆元, 条件是模数必须是质数

    理论支撑:
    设模数为(m), 则对于一个数(i)

    [m = pi +r ag{r<i} ]

    [pi+r equiv 0 (mod m) ]

    [r equiv -pi (mod m) ]

    [1 equiv -pir^{-1} (mod m) ]

    [i^{-1} equiv -pr^{-1} (mod m) ]

    由于(r<i), 故通过这种方法可以线性递推([1,m-1])的逆元。
    (再次强调此方法应用的条件是模数为质数)

    代码():

    int n,p;//n<p, p为模数
    long long inv[maxn];//inv[i]为i在膜p意义下的逆元
    void solve() {
    	inv[1]=1; 
    	for(int i=2;i<=n;++i) inv[i]=(p-p/i)*inv[p%i]%p;
    }
    

    求解同余方程组

    就是求

    [egin{cases} x equiv a_1 (mod m_1)\ x equiv a_2 (mod m_2)\ vdots \ x equiv a_n (mod m_n) end{cases} ]

    的一个解(x)(通常是求最小正整数解)。


    模数两两互质的情况。

    这种情况下一定有解且求解比较简单。

    这种情况下可以为每个(i)构造一个数(node_i), 满足

    [node_i equiv 1(mod m_i) ]

    [node_i equiv 0(mod m_j) ag{$j eq i$} ]

    (x)的一个可行解就是(sum_{i=1}^n a_i *node_i),另外,(x)的通解是

    [sum_{i=1}^n a_i *node_i + k * prod_{i=1}^n m_i ag{$k in Z$} ]

    如何构造(node_i)

    [M_i = frac{prod_{j=1}^n m_j}{m_i} ]

    [I_i equiv M_i^{-1} ag{$mod m_i$} ]

    [node_i equiv M_i * I_i ag{$mod prod_{j=1}^{n} m_j$} ]

    代码():

    #define li long long
    int n;
    li a[maxn], m[maxn];
    //a、m意义见上文
    li CRT() {
    	li M = 1ll, res = 0ll;
    	// res就是要求解的x
    	for(int i=1; i<=n; ++i) M *= m[i]; //这就是M的意义
    	for(int i=1; i<=n; ++i) {
    		li Mi = M/m[i];
    		li Ii = inv(Mi, m[i]);// inv(a, b)是a在膜b意义下的逆元
    		res += Mi * Ii % M * a[i] % M;
    		res %= M;
    	}
    	return ((res%M+M)%M);
    }
    

    模数不一定两两互质的情况

    (X_k)为前(k)个方程组成的方程组的一个特解, 则前(k)个方程组成的方程组的通解为

    [X_k + q * LCM_{i=1}^k m_i ag{$q in Z$} ]

    若前(k+1)个方程组成的方程组有解, 则方程

    [X_k + q * LCM_{i=1}^k m_i equiv a_{k+1} ag{$mod m_{k+1}$} ]

    有解。

    由此, 可以一步一步推出整个方程组的解(有解的话)。

    代码实现中有不少细节, 主要是关于解同余方程时对正负数的处理。

    无解的情况就是exgcd无解的情况。

    代码()(那道板子不需要判断无解):

    // 各个变量基本都和上面的一样
    int n;
    li a[maxn], m[maxn];
    li EXCRT() {
    	li X = a[1], LCM = m[1];
    	for(int i=2; i<=n; ++i) {
    		// k*LCM +m[i]y = c (mod m[i])
    		li c = ((a[i]-X)%m[i] + m[i])%m[i];
    		li x, y, d;
    		exgcd(LCM, m[i], x, y, d);
    		x = slow_mul(x, c/d, m[i]/d); //x 修正为最小正整数解; slow_mul(a,b,c) = (  (a*b)%c +c  )  %  c;
    		X += x*LCM;
    		LCM = LCM/d * m[i];
    		X = (X%LCM+LCM)%LCM;
    	}
    	return X;
    }
    

    计算膜意义下的组合数

    就是计算

    [C_n^m ag{$mod p$} ]


    (p)为质数的情况。
    有个定理叫做(Lucas)定理, 其内容为:

    [C_n^m equiv C_{lfloor n/p floor}^{lfloor m/p floor} * C_{n \% p}^{m \% p} ag{$mod p$,$p$为质数} ]

    代码():

    //jc[i]是预处理的i的阶乘
    li jc[maxp];
    li C(li n, li m, li p) {
    	if(n<m) return 0ll;
    	return jc[n] * ksm(jc[m],p-2,p) % p * ksm(jc[n-m],p-2,p) % p;
    }
    li lucas(li n, li m, li p) {
    	if(!m) return 1ll;
    	if(n<m) return 0ll;
    	return lucas(n/p,m/p,p) * C(n%p,m%p,p) % p;
    }
    

    (p)不一定为质数的情况
    此时可以用(exLucas)求解。

    目标: 求

    [C_n^m mod a ]

    算法过程:

    (a)质因数分解, 即求出(p_1^{c_1}*p_2^{c_2}*cdots*p_k^{c_k} = a)
    分别计算 (ans_j = C_n^m mod p_j^{c_j} (1 leq j leq k))

    最后求出

    [egin{cases} X equiv ans_1 (mod p_1^{c_1}) \ X equiv ans_2 (mod p_2^{c_2}) \ vdots \ X equiv ans_k (mod p_k^{c_k}) end{cases} ]

    的解(X)(X)即为(C_n^m mod a)

    如何求解(C_n^m mod p^k (p为质数))

    首先, (C_n^m = frac{n!}{(n-m)!m!}), 但((n-m)!m!)中可能有(p)这个因子, 此时((n-m)!m!)在膜(p)意义下没有逆元, 于是就要将((n-m)!m!)中的(p)因子提出来, 最后再乘回去。 但是提出来的(p)因子是(frac{1}{p})形式的, 还是无法计算, 此时就要将(n!)(p)因子也提出来。

    (n!)含有(p)因子的个数为(num_n)(m!)含有(p)因子的个数为(num_m)((n-m)!)含有(p)因子的个数为(num_{n-m})

    似乎可以证明(p^{num_n-num_m-num_{n-m}} geq 0)。(我不会证也不想证qwq)

    于是最后要算的就是

    [frac{frac{n!}{p^{num_n}}}{frac{m!}{p^{num_m}} frac{(n-m)!}{p^{num_{n-m}}}} * p^{num_n-num_m-num_{n-m}} mod p^k ]

    如何求解(frac{n!}{p^{num_n}} mod p^k)

    首先将(n!)展开为

    [1 * 2 *3 * 4 * cdots * n-1 * n ]

    将其中(p)的倍数都提出来

    [(p * 2p *3p * 4p * cdots * lfloor n/p floor p) * 1 * cdots * unknown ]

    再操作一下

    [Big( p^{lfloor n/p floor} Big) * Big( 1 * 2 *3 * 4 * cdots * lfloor n/p floor Big) * Big( 1 * cdots * unknown(都是与p互质的数) Big) ]

    将第一部分扔掉(因为要求的是(frac{n!}{p^{num_n}} mod p^k)),至此可以发现问题被缩小了。

    但是如何求([1,n])中与(p)互质的数的乘积(mod p^k)的值呢?

    首先(x equiv x + p^k (mod p^k)), 然后就可以(O(p^k))([1,n])中与(p)互质的数的乘积(mod p^k)的值了。

    代码()(完整)(原题链接):

    #include<bits/stdc++.h>
    using namespace std;
    #define li long long
    li ksm(li a, li b, li p) {
    	li res = 1ll;
    	for(;b;b>>=1, a=(a*a)%p)
    		if(b&1) res = (res*a)%p;
    	return res%p;
    }
    void exgcd(li a, li b, li &x, li &y) {
    	!b ? x=1,y=0 : (exgcd(b,a%b,y,x), y-=x*(a/b));
    }
    li inv(li a, li p) {
    	// ax + py = 1;
    	li x, y;
    	exgcd(a, p, x, y);
    	x = ((x%p+p)%p);
    	return x;
    }
    
    li fac(li n, li p, li pk) {
    	if(!n) return 1ll;
    	li res = 1ll;
    	for(li i=0; i<pk; ++i) if(i%p) res = res*i % pk;
    	res = ksm(res, n/pk, pk);
    	for(li i=0; i<=n%pk; ++i) if(i%p) res = res*i % pk;
    	return res * fac(n/p,p,pk) % pk;
    }
    
    li C(li n, li m, li p, li pk) {
    	if(!m) return 0ll;
    	li f1 = fac(n,p,pk), f2=fac(m,p,pk), f3=fac(n-m,p,pk), cnt = 0ll;
    	for(li i=n;i;i/=p) cnt += i/p;
    	for(li i=m;i;i/=p) cnt -= i/p;
    	for(li i=n-m;i;i/=p) cnt -= i/p;
    	return f1 * inv(f2,pk) % pk * inv(f3,pk) % pk * ksm(p,cnt,pk) % pk;
    }
    
    int tot;
    li ri[10001], riri[10001];
    li CRT() {
    	li M = 1ll, res = 0ll;
    	for(int i=1; i<=tot; ++i) M *= riri[i];
    	for(int i=1; i<=tot; ++i) {
    		int nd = M/riri[i];
    		res += ri[i] * nd % M * inv(nd,riri[i]) % M;
    		res %= M;
    	}
    	return ((res%M+M)%M);
    }
    li exlucas(li n, li m, li p) {
    	for(li i=2; i<=sqrt(p); ++i) {
    		li md = 1ll;
    		while(p%i==0) p/=i, md*=i;
    		if(md>1) ri[++tot] = C(n,m,i,md), riri[tot] = md;
    	}
    	if(p>1) ri[++tot] = C(n,m,p,p), riri[tot] = p;
    	return CRT();
    }
    int main() {
    	li n, m, p;
    	cin >> n >> m >> p;
    	cout << exlucas(n,m,p);
    	return 0;
    }
    

    原根!

    从别人课件扒下来的

    定义 (g) 是质数 (p) 的原根, 当且仅当 (1,g,g^2, cdots,g^{p-2}) 互不相同 (当然是膜 (p) 意义下的)。

    正整数其实也有原根来着

    定理 : 任何质数都有原根

    原根的个数不一定是唯一的。

    如何判定一个数 (g) 是否是 (p) 的原根?
    (g) 不是 (p) 的原根, 则存在一正整数 (0 le k < p-1) 使得
    (g^k equiv 1 (mod p))
    如此,
    再由 (gcd(k, p-1) = xk + y(p-1))
    得出

    [g^{gcd(k,p-1)} equiv 1 ag{$mod p$} ]

    这样, 就可以枚举 (p-1) 的真因子来判断 (g) 是否为 (p) 的原根。
    设存在一个 (p-1) 的真因子 (s) 使得 (g^s equiv 1(mod p)), 则 (g) 不是 (p) 的原根, 反之 (g)(p) 的原根。
    其实只要判断形如 "((p-1)/一个质数, 这个质数|(p-1))" 形式的因子就好了。

    如何找质数 (p) 的原根?
    2 开始枚举啊
    (2) 枚举即可, 一般来说找的挺快,不会解释。

    代码?

    //找质数 p 的原根 g 
    //和原根表对了几个,应该阔以保证正确,但是 p 太大会炸, 不知道为啥 
    
    #include<bits/stdc++.h>
    using namespace std;
    #define li long long
    li ksm(li a, li b, li p) {
    	li res = 1ll;
    	for(;b;b>>=1, a=(a*a)%p)
    		if(b&1) res=(res*a)%p;
    	return res%p;
    }
    const int maxn = 1000005;
    
    li p, g;
    li tot, P[maxn], C[maxn];
    li factot, fac[maxn];
    
    void getfac(li now, li nownum) {
    	if(now==tot+1) {
    		fac[++factot] = nownum;
    		return;
    	}
    	for(int i=0; i<=C[now]; ++i) {
    		getfac(now+1, nownum);
    		nownum *= P[now];
    	}
    }
    bool jd() {
    	for(int i=1;i<=factot;++i) if(fac[i]!=p-1 && ksm(g,fac[i],p) == 1ll) {
    		return false;
    	}
    	return true;
    }
    
    void gtfac(li x) {
    	for(li i=2; i<=sqrt(x); ++i) if(x%i==0) {
    		P[++tot] = i, C[tot] = 0;
    		while(x%i==0) ++C[tot], x/=i;
    	}
    	if(x>1) P[++tot] = x, C[tot]=1;
    	getfac(1,1ll);
    }
    int main()
    {
    	cin >> p;
    	gtfac(p-1);
    //	for(int i=1;i<=factot;++i) cout<<fac[i]<<'
    ';
    //	g=3;
    //	cout<<jd();
    	cout<<p<<'
    ';
    	if(p==3) {
    		cout<<2;
    		return 0;
    	}
    	for(g=2ll; g<p-1; g+=1ll) if(jd()){
    		cout << g;
    		break;
    	}
    	return 0;
    }
    

    [gcd(x_1, cdots, x_n) Rightarrow O(n + log w) ]

    证明:

    每做两次辗转相除, gcd/2, 这样的次数最多只有 log w 次。


    小凯的疑惑

    ax + by = d

    由于 gcd(a,b) = 1, 故这个方程一定有整数解

    假设 d 已经固定, 则一组解 (x,y) 可以看成平面直角系上的一个点。

    根据解 (x,y), 得出通解 (x+kb, y-ka), (k 是整数), 把这个解集里的点都放在直角坐标系上, 它们在一条斜率为负数且斜率一定的直线上, 直线上相邻两个点的坐标差完全一样;且既然斜率是负数, 那么这条直线一定穿过 2、4 象限, 而且由于 d 是正数, 这条直线一定过第一象限。

    其实通过前面的分析, 可以看出直线与 y 的交点的纵坐标 h 无限增大的时候, 一定存在一个 h, 使得直线上的解一定在第一象限或 x 轴正半轴或 y 正半轴上, 那么就猜测 d 增大的时候, 直线与 y 轴的交点的纵坐标增大

    还是考虑那个 横坐标+b, 纵坐标-a 就从第二象限跳到第四象限的点,由于 d 要尽量大, 所以 x = -1, 跳之后 x = b-1,无法继续考虑。

    由于 d 要尽量大, 所以令

  • 相关阅读:
    Run Book Automation
    Android.mk中的经常使用语法
    层的匀速运动波动处理
    【ArcGIS 10.2新特性】ArcGIS 10.2 for Server常见问题
    WPF-19:分享一个样式(左右滑动选中的checbox)
    [置顶] 程序员期望月薪那些事儿
    【VB/.NET】Converting VB6 to VB.NET 【Part II】【之四】
    两种方式给列表增加自动增长序号列
    在后台运行erlang;在需要时连回交互模式
    php设计模式——UML类图
  • 原文地址:https://www.cnblogs.com/tztqwq/p/12735627.html
Copyright © 2011-2022 走看看