dalao博客,至少很好看。。
因为本人数论实在渣渣,但是考试确是得考的,只好尽早学,尽早掌握。
最大公因数
普通gcd
O(log(min(a,b)))
1 inline int gcd(int x,int y) 2 { 3 return y == 0 ? x : gcd(y, x % y) 4 }
二进制优化gcd
O(log(min(a,b))) 没有%常数减小
1 inline int bsgcd(int x, int y) 2 { 3 if(x == y) return x; 4 if(x < y) x ^= y ^= x ^= y; 5 if(!(x & 1)) //x偶 y偶 gcd(x,y)=2*gcd(x/2,y/2),x偶 y奇 gcd(x,y)=gcd(x/2,y) 6 return (!(y & 1)) ? bsgcd(x >> 1, y >> 1) << 1 : bsgcd(x >> 1, y); 7 //x奇 y偶 gcd(x, y)=gcd(x,y/2) 8 return (!(y & 1)) ? bsgcd(x, y >> 1) : bsgcd(y, x - y); 9 }
素数判定
埃拉托色尼筛法
O(nloglogn)
不够优越,所以不学。
线性筛
O(n)
1 const int MAXN = 10000001; 2 3 bool notpri[MAXN]; 4 int n, m, cnt, prime[MAXN]; 5 6 inline void get_prime() 7 { 8 int i, j; 9 notpri[0] = notpri[1] = 1; // 1 表示不是质数 10 for(i = 2; i <= n; i++) 11 { 12 if(!notpri[i]) prime[++cnt] = i; 13 for(j = 1; i * prime[j] <= n; j++) 14 { 15 notpri[i * prime[j]] = 1; 16 if(!(i % prime[j])) break; 17 } 18 } 19 }
Miller-Rabin 大素数判定
http://www.cnblogs.com/vongang/archive/2012/03/15/2398626.html
O(log32n)
费马小定理:对于素数p和任意整数a,有ap ≡ a(mod p)(同余)。反过来,满足ap ≡ a(mod p),p也几乎一定是素数。
伪素数:如果p是一个正整数,如果存在和p互素的正整数a满足 ap-1 ≡ 1(mod p),我们说p是基于a的伪素数。如果一个数是伪素数,那么它几乎肯定是素数。
Miller-Rabin测试:不断选取不超过p-1的基a(s次),计算是否每次都有ap-1 ≡ 1(mod p),若每次都成立则p是素数,否则为合数。
还有一个定理,能提高Miller测试的效率:
二次探测定理:如果p是奇素数,则 x2 ≡ 1(mod p)的解为 x = 1 || x = p - 1(mod p);
这个算法的正确率据说是3/4,所以多测几次。。
1 #include <cstdio> 2 #include <ctime> 3 #include <cstdlib> 4 #define LL long long 5 6 int n, ans; 7 const int S = 5; 8 9 //快速计算 a * b % p 10 inline LL mod_mul(LL a, LL b, LL p) 11 { 12 LL ret = 0; 13 while(b) 14 { 15 if(b & 1) ret = (ret + a) % p; 16 a = (a << 1) % p; 17 b >>= 1; 18 } 19 return ret; 20 } 21 22 //快速计算 a ^ b % p 23 inline LL mod_exp(LL a, LL b, LL p) 24 { 25 LL ret = 1; 26 while(b) 27 { 28 if(b & 1) ret = mod_mul(ret, a, p); 29 a = mod_mul(a, a, p); 30 b >>= 1; 31 } 32 return ret; 33 } 34 35 inline bool miller_rabin(LL p) 36 { 37 if(p == 2 || p == 3 || p == 5 || p == 7 || p == 11) return 1; 38 if(p == 1 || !(p % 2) || !(p % 3) || !(p % 5) || !(p % 7) || !(p % 11)) return 0; 39 LL x, pre, u; 40 int i, j, k = 0; 41 u = p - 1; // 要求 u = p - 1 , (x ^ u) % p 42 while(!(u & 1)) k++, u >>= 1; // 如果 u 为偶数则右移,k 记录位数 43 srand((LL)time(0)); 44 for(i = 0; i < S; i++) 45 { 46 x = rand() % (p - 2) + 2; // 在 [2, p) 中选取随机数 47 if(!(x % p)) continue; // 如果不是互质 48 x = mod_exp(x, u, p); // 先计算 (x ^ u) % p 49 pre = x; 50 for(j = 0; j < k; j++) // 把移位补 51 { 52 x = mod_mul(x, x, p); // 运用二次探测 53 if(x == 1 && pre != 1 && pre != p - 1) return 0; 54 // 如果 (x ^ 2) % p == 1 但是 x != 1, x != p - 1, 就不是素数 55 pre = x; 56 } 57 if(x != 1) return 0; // 费马小定理 58 } 59 return 1; 60 } 61 62 int main() 63 { 64 LL x; 65 while(~scanf("%d", &n)) 66 { 67 ans = 0; 68 while(n--) 69 { 70 scanf("%lld", &x); 71 if(miller_rabin(x)) ans++; 72 } 73 printf("%d ", ans); 74 } 75 return 0; 76 }
拓展欧几里得
能求出 ax + by = gcd(a,b) 的一组解,还能顺带求出 gcd
正常版
1 inline int exgcd(int a, int b, int &x, int &y) { 2 if (!b) {x = 1, y = 0; return a;} 3 int r = exgcd(b, a % b, x, y); 4 int t = x, x = y, y = t - a / b * y; 5 return r; 6 }
极简版
1 inline int exgcd(int a, int b, int &x, int &y) 2 { 3 if(!b){x = 1, y = 0; return a;} 4 int r = exgcd(b, a % b, y, x); 5 y -= a / b * x; 6 return r; 7 }
证明——其实就是推一下式子,
ax+by=gcd(a,b) ①
bx'+(a%b)y'=gcd(b,a%b) ②
已知x'和y'求x和y
∵gcd(a,b)=gcd(b,a%b)
∴ax+by=bx'+(a%b)y'
ax+by=bx'+(a-a/b*b)y'
ax+by=bx'+ay'-a/b*by'
ax+by=ay'+b(x'-a/b*y')
∴x=y',y=x'-a/b*y'
证毕,但要注意除法都是下取整的
快速幂
递归版
1 #define LL long long 2 inline LL pow(LL a, LL b, LL p) 3 { 4 if(!b) return 1; 5 LL ans = pow(a, b >> 1, p); 6 ans = ans * ans % p; 7 if(b & 1) ans = a * ans % p; 8 return ans; 9 }
非递归版
1 #define LL long long 2 inline LL pow(LL a, LL b, LL p) 3 { 4 LL ans = 1; 5 while(b) 6 { 7 if(b & 1) ans = ans * a % p; 8 a = a * a % p; 9 b >>= 1; 10 } 11 return ans; 12 }
膜意义下的逆元
若 a*x≡1(mod b),且a和b互质,则称x为a的逆元,记作a-1
逆元的求解方法:
1.快速幂
如果p是质数,那么 a 的逆元为 ap-2,好像是通过什么费马小定理证明的,记住就行
2.扩展欧几里得
定义式可以转化为a*x+b*y==1,所以自然可以利用扩展欧几里得求啦!
3.线性递推
首先1-1≡1(mod p),
假设p=k*i+r,r<i,1<i<p,那么k*i+r≡0(mod p)
等式两边同时乘上i-1,r-1,式子就变成k*r-1+i-1≡0(mod p)
i-1≡-k*r-1
i-1≡-(p/i)*(p mod i)-1
这样便可以线性递推求解
#include <cstdio> #define N 3000001 int n, p; int inv[N] = {0, 1}; int main() { int i; scanf("%d %d", &n, &p); puts("1"); for(i = 2; i <= n; i++) { inv[i] = 1ll * -(p / i) * inv[p % i] % p; printf("%d ", (inv[i] + p) % p); } return 0; }
卢卡斯定理——用于解决组合数的取膜问题
卢卡斯定理
C(n,m)%p=C(n%p,m%p)*C(n/p,m/p)%p
至于证明。。算了吧,背过得了
因为p不大
所以C(n%p,m%p)可直接求出,
但是如果n%p<m%p时返回0
C(n/p,m/p)可继续递归求解
当m==0时返回1
#include <cstdio> #define N 200001 #define LL long long int T; int F[N], inv[N]; inline int C(int n, int m, int p) { if(m > n) return 0; return (LL)F[n] * inv[m] % p * inv[n - m] % p; } inline int Lucas(int n, int m, int p) { if(!m) return 1; return (LL)C(n % p, m % p, p) * Lucas(n / p, m / p, p) % p; } int main() { int i, n, m, p; scanf("%d", &T); while(T--) { scanf("%d %d %d", &n, &m, &p); n += m; F[0] = F[1] = inv[0] = inv[1] = 1; for(i = 2; i <= p; i++) inv[i] = -(LL)(p / i) * inv[p % i] % p; for(i = 2; i <= p; i++) F[i] = (LL)F[i - 1] * i % p, inv[i] = (LL)inv[i] * inv[i - 1] % p; printf("%d ", (Lucas(n, m, p) + p) % p); } return 0; }