\[TOC\]
# 数论模板
## 一, 线性筛
所有数组定义在外面, 自动置为全0即可
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 void init () { cnt = 0 ; is_prime[1 ] = 1 ; for (int i = 2 ; i <= maxn; i++) { if (!is_prime[i]) { prime[++cnt] = i; } for (int j = 1 ; j <= cnt; j++) { if (prime[j] * i >= maxn) break ; is_prime[i * prime[j]] = 1 ; if (i % prime[j] == 0 ) break ; } } }
### 1.分解质因数
通过素数筛法打表, 然后只能通过优化暴力, **最快了**
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 int get_primefactors (int n) { int sum = 0 ; for (int i = 1 ; i <= tot && pri[i] * pri[i] <= n; ++i) { if (n % pri[i] == 0 ) { fac[++sum] = pri[i]; while (n % pri[i] == 0 ) { n /= pri[i]; } } } if (n > 1 ) { fac[++sum] = n; } return sum; }
### 2.区间线性筛
a, b可以很大, 但是b - a 的结果可以接受
模板 : 输出a 到 b之间的素数个数
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 #define ll long long const int MAXN = 100005 ;bool is[MAXN];int pri[MAXN / 10 ];int cnt;bool ans[MAXN];void init () { memset (is, true , sizeof (is)); is[0 ] = is[1 ] = false ; cnt = 0 ; for (int i = 2 ; i < MAXN; i++) { if (is[i]) { pri[++cnt] = i; } for (int j = 1 ; j < cnt && pri[j] * i < MAXN; j++) { is[pri[j] * i] = false ; if (i % pri[j] == 0 ) { break ; } } } } int solve (ll a, ll b) { memset (ans, true , sizeof (ans)); for (int i = 1 ; i <= cnt && 1ll * pri[i] * pri[i] <= b; i++) { ll s = (a + pri[i] - 1 ) / pri[i]; if (s < 2 ) { s = 2 ; } s *= pri[i]; for (; s <= b; s += pri[i]) { ans[s - a] = false ; } } int res = 0 ; for (int i = 0 ; i <= b - a; i++) { if (ans[i]) { res++; } } if (a == 1 ) { res--; } return res; } int main () { init (); int a, b; cin >> a >> b; cout << solve (a, b); }
### 3.求各种积性函数
YC究极版!
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 #include <bits/stdc++.h> #define ll long long using namespace std;const int N = 1e7 + 5 ;bool is[N]; int pri[N]; int phi[N]; int mu[N]; int e[N]; int d[N]; int fac[N]; int tot; void init () { memset (is, true , sizeof (is)); tot = 0 ; is[0 ] = is[1 ] = 0 ; mu[1 ] = 1 ; phi[1 ] = 1 ; fac[1 ] = 1 ; e[1 ] = 0 ; d[1 ] = 1 ; for (int i = 2 ; i < N; i++) { if (is[i]) { pri[++tot] = i; phi[i] = i - 1 ; mu[i] = -1 ; e[i] = 1 ; d[i] = 2 ; fac[i] = i; } for (int j = 1 ; j <= tot && pri[j] * i < N; j++) { is[pri[j] * i] = 0 ; if (i % pri[j] == 0 ) { phi[i * pri[j]] = phi[i] * pri[j]; mu[i * pri[j]] = 0 ; d[pri[j] * i] = d[i] / (e[i] + 1 ) * (e[i] + 2 ); e[pri[j] * i] = e[i] + 1 ; fac[i * pri[j]] = fac[i] * pri[j]; break ; } else { mu[i * pri[j]] = -mu[i]; phi[i * pri[j]] = phi[i] * (pri[j] - 1 ); e[pri[j] * i] = 1 ; d[pri[j] * i] = d[i] * d[pri[j]]; fac[i * pri[j]] = pri[j]; } } } } int main () { init (); int a, b; cin >> a >> b; cout << d[a] << ' ' << d[b]; return 0 ; }
## 二, Min\_25筛
#### 求十亿(1e9)内所有质数的和
$O(n ^{3/4})$
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 const int N = 1000010 ;typedef long long LL;namespace Min25 { int prime[N], id1[N], id2[N], flag[N], ncnt, m; LL g[N], sum[N], a[N], T, n; inline int ID (LL x) { return x <= T ? id1[x] : id2[n / x]; } inline LL calc (LL x) { return x * (x + 1 ) / 2 - 1 ; } inline LL f (LL x) { return x; } inline void init () { T = sqrt (n + 0.5 ); for (int i = 2 ; i <= T; i++) { if (!flag[i]) prime[++ncnt] = i, sum[ncnt] = sum[ncnt - 1 ] + i; for (int j = 1 ; j <= ncnt && i * prime[j] <= T; j++) { flag[i * prime[j]] = 1 ; if (i % prime[j] == 0 ) break ; } } for (LL l = 1 ; l <= n; l = n / (n / l) + 1 ) { a[++m] = n / l; if (a[m] <= T) id1[a[m]] = m; else id2[n / a[m]] = m; g[m] = calc (a[m]); } for (int i = 1 ; i <= ncnt; i++) for (int j = 1 ; j <= m && (LL)prime[i] * prime[i] <= a[j]; j++) g[j] = g[j] - (LL)prime[i] * (g[ID (a[j] / prime[i])] - sum[i - 1 ]); } inline LL solve (LL x) { if (x <= 1 ) return x; return n = x, init (), g[ID (n)]; } } int main () { LL n; scanf ("%lld" , &n); printf ("%lld\n" , Min25::solve (n)); }
#### 求积性函数的前缀和
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 #include<bits/stdc++.h> using namespace std; #define ll long long const int maxn = 2000000+100; /****************************** f()函数中(31-37行) 填函数在质数幂次处的表达式 pow_sum()函数中(38-43行) 填幂和函数(如果需要更高次的话可以在这里添加) 202-205行按要求填写 f_p[][0/1/2/3/...]分别代表质数个数/质数和/质数平方和/质数三次方和/...根据自己需要添加 例:如果该函数在质数处表达式为f(p) = p^2+3*p+1 则表明需要质数个数/质数和/质数平方和 即f_p[][0],f_p[][1],f_p[][2] ******************************/ ll poww(ll a,ll b){ ll res = 1; ll base = a; while(b){ if(b&1){ res *= base; //res %= mod; } base *= base; //base %= mod; b>>=1; } return res; } inline ll f(ll p,int e){ if(p==1||e==0) return 1; ///return f(p^e) ll res = poww(p,e); return res*res+3*res+1; } ll pow_sum(ll n,int k){ ///return sum(i^k),i from 1 to n. if(k==0) return n; if(k==1) return n*(n+1)/2; if(k==2) return n*(n+1)*(2*n+1)/6; } ll f_p[maxn][3];///F_prime(id(n/i)) ll n; int n_2; ///(int)sqrt(n) int n_3; ///(int)pow(n,1.0/3.0) int n_6; ///(int)pow(n,1.0/6.0) ll val_id[maxn]; ///give the id, return the id-th number like 'n/i' ,(val_id[1] = 1) int val_id_num; ///how many numbers like 'n/i' int val_id_num_3; ///how many numbers like 'n/i' below n/n_3; int p[200000+100]; bool isp[maxn]; int p_sz_2; ///pi(n_2) int p_sz_3; ///pi(n_3) int p_sz_6; ///pi(n_6) void init(){ n_2 = (int)sqrt(n); n_3 = (int)pow(n,1.0/3.0); n_6 = (int)pow(n,1.0/6.0); val_id_num = 0; for(ll i=1;i<=n;){ val_id[++val_id_num] = i; if(i==n) break; i = n/(n/(i+1)); } memset(isp,1,sizeof isp); isp[1] = 0; for(int i=2;i<=n_2;i++){ if(isp[i]){ p[++p_sz_2] = i; if(i<=n_3) p_sz_3++; if(i<=n_6) p_sz_6++; } for(int j=1;j<=p_sz_2&&p[j]*i<=n_2;j++){ isp[i*p[j]] = 0; if(i%p[j]==0) break; } } } inline int get_id(ll k){ ///give a number like 'n/i', return the id of it if(k>n_2) return val_id_num-n/k+1; else return k; } ll c[maxn]; int lowbit(int n){return n & (-n);} void add(int x,ll d){ while(x<maxn){ c[x]+=d; x+=lowbit(x); } } ll sum(int x){ ll ans=0; while(x){ ans+=c[x]; x-=lowbit(x); } return ans; } struct node{ int k_max; ll val; ll f_val; }; void update_bfs(int k,int type){ queue<node> q; while(!q.empty()) q.pop(); int e = 1; for(ll i=p[k];i<n/n_3;i*=p[k]){ node st; st.k_max = k; st.val = i; if(type==-1)st.f_val = f(p[k],e); else st.f_val = poww(i,type); q.push(st); e++; } while(!q.empty()){ node hd = q.front(); q.pop(); if((hd.val!=p[hd.k_max]&&type>=0)||type==-1) {//if(type==-1)cout << "****" << hd.val << "****" << hd.f_val << endl; ll w = n/hd.val; w = n/w;//cout << hd.val << "[" << w<<" , " << val_id[val_id_num] << "]" << endl; if(type==-1){ add(get_id(w),hd.f_val); add(val_id_num+1,-1ll*hd.f_val); } else{ add(get_id(w),-1ll*hd.f_val); add(val_id_num+1,hd.f_val); } } for(int i=hd.k_max+1;hd.val*p[i]<n/n_3&&i<=p_sz_2;i++){ ll res = p[i]; for(int e=1;;e++){ if(hd.val*res<n/n_3){ node nxt; nxt.k_max = i; nxt.val = hd.val*res; if(type==-1) nxt.f_val = hd.f_val*f(p[i],e); else nxt.f_val = hd.f_val*poww(res,type); q.push(nxt); } else break; res *= p[i]; } } } } void get_f_p(ll n,int times){ for(int i=1;i<=val_id_num;i++){ for(int j=0;j<=times;j++){ f_p[i][j] = pow_sum(val_id[i],j)-1; } } int now; //for(now=1;now<=p_sz_2;now++){ for(now=1;p[now]<=n_6;now++){ for(int j=val_id_num;j>=1;j--){ ll w = val_id[j]/p[now]; if(w<p[now]) break; ll val=1; for(int k = 0;k<=times;k++){ f_p[j][k] = f_p[j][k] - val*(f_p[get_id(w)][k]-f_p[p[now-1]][k]); val *= p[now]; } } } int nnow = now; int val = 1; for(int tt = 0;tt<=times;tt++){ now = nnow; memset(c,0,sizeof c); add(1,f_p[1][tt]); for(int i=2;val_id[i]<n/n_3;i++){ add(i,f_p[i][tt] - f_p[i-1][tt]); } for(;p[now]<=n_3;now++){ for(int j=val_id_num;j>=1;j--){ ll w = val_id[j]/p[now]; if(val_id[j]<n/n_3) break; if(w<p[now]) break; if(w<n/n_3) f_p[j][tt] = f_p[j][tt] - (sum(get_id(w)) - sum(p[now-1]))*poww(p[now],tt); else f_p[j][tt] = f_p[j][tt] - (f_p[get_id(w)][tt]-sum(p[now-1]))*poww(p[now],tt); } update_bfs(now,tt); } for(int i=1;i<=val_id_num&&val_id[i]<n/n_3;i++) f_p[i][tt] = sum(i); for(;now<=p_sz_2;now++){ for(int j=val_id_num;j>=1;j--){ ll w = val_id[j]/p[now]; if(val_id[j]<n/n_3) break; if(w<p[now]) break; f_p[j][tt] -= (f_p[get_id(w)][tt]-f_p[p[now-1]][tt])*poww(p[now],tt); } } } for(int i=1;i<=val_id_num;i++){ ///if f(p) = p^2+3p+1,then write:f_p[i][0] = f_p[i][2] + 3*f_p[i][1] + f_p[i][0]; f_p[i][0] = f_p[i][2] + 3*f_p[i][1] + f_p[i][0]; } } ll F[2000000+100]; void get_f_3(ll n){ ///V(F_{pi(n^(1/3))+1},n) ll q = p[p_sz_3+1]; for(int now=1;now<=val_id_num;now++){ if(val_id[now]<q){ F[now] = 1; } else if(val_id[now]<q*q){ F[now] = 1+(f_p[now][0]-f_p[q-1][0]); } else{ F[now] = 1+(f_p[now][0]-f_p[q-1][0]); for(int pp=p_sz_3+1;p[pp]<=(int)(sqrt(val_id[now]))&&pp<=p_sz_2;pp++){ F[now] += f(p[pp],2) + (f(p[pp],1))*(f_p[get_id(val_id[now]/p[pp])][0]-f_p[get_id(p[pp])][0]); } } } } void get_f_6(ll n){ ///V(F_{pi(n^(1/6))+1},n) memset(c,0,sizeof c); add(1,F[1]); for(int i=2;val_id[i]<n/n_3;i++){ add(i,F[i] - F[i-1]); } for(int k=p_sz_3;k>p_sz_6;k--){ int now = val_id_num; for(;val_id[now]>=n/n_3;now--){ int e = 1; ll _p = p[k]; while(val_id[now]/_p){ if(val_id[now]/_p>=n/n_3){ F[now] += F[get_id(val_id[now]/_p)]*f(p[k],e); } else{ F[now] += sum(get_id(val_id[now]/_p))*f(p[k],e); } _p *= p[k]; e++; } } if(k==1) break; //cout << "******" << p[k] << "******" << n/n_3 << endl; update_bfs(k,-1);///bfs to update [lpf(i)==P{k-1}]f(i) } for(int i=1;i<=val_id_num&&val_id[i]<n/n_3;i++) F[i] = sum(i); } void get_f(ll n){ for(int k=p_sz_6;k>=1;k--){ for(int now = val_id_num;now>=1;now--){ int e = 1; ll _p = p[k]; while(val_id[now]/_p){ F[now] += F[get_id(val_id[now]/_p)]*f(p[k],e); _p *= p[k]; e++; } } } } int main(){//n = 1000000000; //1e10:455052511,0.83s/0.58s 1e12:37607912018 9.224s/5.105s cin >> n; init(); get_f_p(n,2); get_f_3(n); get_f_6(n); get_f(n); for(int i=1;i<=val_id_num;i++){ cout << val_id[i] << " : " << F[i] << endl; } } ``` ## 三, 唯一分解定理(算术基本定理) 暴力分解即可 代码没有 ## 四, GCD/EXGCD 用途: **ax + by = gcd(a,b) 方程求解** ```C++ #include <bits/stdc++.h> using namespace std; int a, b, x, y; void EXGCD(int a, int b, int &x, int &y) { if (!b) { x = 1; y = 0; return; } EXGCD(b, a%b, y, x); y -= a / b * x; } int main() { scanf("%d %d", &a, &b); EXGCD(a,b,x,y); printf("%d %d", x, y); return 0; }
使用STL的代码:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 #include <bits/stdc++.h> using namespace std;#define pir pair<int,int> int a, b, x, y;pir EXGCD (int a,int b) { if (!b) return make_pair (1 , 0 ); pir tmp = EXGCD (b, a%b); return make_pair (tmp.second ,tmp.first - a/b*tmp.second); } int main () { scanf ("%d %d" , &a, &b); pir ans = EXGCD (a,b); printf ("%d %d" ,ans.first,ans.second); return 0 ; }
## 五, 中国剩余定理(孙子定理)
### 1.CRT(保证$m\_1,m\_2…m\_n$之间两两互质,求最小的_x_)
用途: **同余方程组求解!**

我们求出了任意一个满足条件的_x_之后,只需要对其模_M_就是最终的答案(因为_M_M是所有数的_lcm_)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 #include <cstdio> #define ll long long void gcd (ll a,ll b,ll &d,ll &x,ll &y) { if (b==0 ){ d=a; x=1 ,y=0 ; } else { gcd (b,a%b,d,y,x); y-=(a/b)*x; } } ll China (int n,ll *m,ll *a) { ll M=1 ,d,y,x=0 ; for (int i=0 ;i<n;i++) M*=m[i]; for (int i=0 ;i<n;i++){ ll w=M/m[i]; gcd (m[i],w,d,d,y); x=(x+y*w*a[i])%M; } return (x+M)%M; }
### 2\. EXCRT($m\_1,m\_2…m\_n$之间不再限制互质,求最小的_x_)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 #include <cstdio> #include <algorithm> #include <iostream> using namespace std;typedef long long ll;const int N=100000 +10 ;int n;ll a[N],r[N]; ll exgcd (ll a,ll b,ll &x,ll &y) { if (b==0 ){ x=1 ,y=0 ;return a; } ll ret=exgcd (b,a%b,y,x); y-=(a/b)*x; return ret; } ll excrt () { ll M=a[1 ],R=r[1 ],x,y,d; for (int i=2 ;i<=n;i++){ d=exgcd (M,a[i],x,y); if ((R-r[i])%d) return -1 ; x=(R-r[i])/d*x%a[i]; R-=M*x; M=M/d*a[i]; R%=M; } return (R%M+M)%M; } int main () { while (scanf ("%d" ,&n)!=EOF){ for (int i=1 ;i<=n;i++){ scanf ("%lld%lld" ,&a[i],&r[i]); } printf ("%lld\n" ,excrt ()); } return 0 ; }
## 六, 米勒罗宾\_大素数测试
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 #include <bits/stdc++.h> #define ll long long using namespace std;const int S = 8 ;ll mult_mod (ll a, ll b, ll c) { a %= c; b %= c; ll ans = 0 ; ll tmp = a; while (b) { if (b & 1 ) { ans += tmp; if (ans > c) { ans -= c; } } tmp <<= 1 ; if (tmp > c) { tmp -= c; } b >>= 1 ; } return ans; } ll pow_mod (ll a, ll n, ll mod) { ll ans = 1 ; ll tmp = a % mod; while (n) { if (n & 1 ) { ans = mult_mod (ans, tmp, mod); } tmp = mult_mod (tmp, tmp, mod); n >>= 1 ; } return ans; } bool check (ll a, ll n, ll x, ll t) { ll ans = pow_mod (a, x, n); ll last = ans; for (int i = 1 ; i <= t; ++i) { ans = mult_mod (ans, ans, n); if (ans == 1 && last != 1 && last != n - 1 ) { return true ; } last = ans; } if (ans != 1 ) { return true ; } return false ; } bool Miller_Rabin (ll n) { if (n < 2 ) { return false ; } if (n == 2 ) { return true ; } if ((n & 1 ) == 0 ) { return false ; } ll x = n - 1 ; ll t = 0 ; t = (ll)log2 (x & (-x)); x = x / (x & (-x)); srand (time (NULL )); for (int i = 0 ; i < S; ++i) { ll a = rand () % (n - 1 ) + 1 ; if (check (a, n, x, t)) { return false ; } } return true ; } int main () { int t, ans; ll tem; scanf ("%d" , &t); ans = 0 ; while (t--) { scanf ("%lld" , &tem); if (Miller_Rabin (tem)) { ans++; } } cout << ans << endl; return 0 ; }
### 增: Pallord\_rho 大整数分解质因数
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 #include <bits/stdc++.h> #define ll long long #define INF 0x3f3f3f3f #define maxn 10000+10 #define cle(a) memset(a, 0, sizeof(a)) const double eps = 1e-5 ;using namespace std;const int S = 30 ;ll mult_mod (ll a, ll b, ll c) { a %= c; b %= c; ll ret = 0 ; while (b) { if (b & 1 ) { ret += a; ret %= c; } a <<= 1 ; if (a >= c)a %= c; b >>= 1 ; } return ret; } ll pow_mod (ll x, ll n, ll mod) { if (n == 1 )return x % mod; x %= mod; ll tmp = x; ll ret = 1 ; while (n) { if (n & 1 ) ret = mult_mod (ret, tmp, mod); tmp = mult_mod (tmp, tmp, mod); n >>= 1 ; } return ret; } bool check (ll a, ll n, ll x, ll t) { ll ret = pow_mod (a, x, n); ll last = ret; for (int i = 1 ; i <= t; i++) { ret = mult_mod (ret, ret, n); if (ret == 1 && last != 1 && last != n - 1 ) return true ; last = ret; } if (ret != 1 ) return true ; return false ; } bool Miller_Rabin (ll n) { if (n < 2 )return false ; if (n == 2 )return true ; if ((n & 1 ) == 0 ) return false ; ll x = n - 1 ; ll t = 0 ; while ((x & 1 ) == 0 ) { x >>= 1 ; t++; } for (int i = 0 ; i < S; i++) { ll a = rand () % (n - 1 ) + 1 ; if (check (a, n, x, t)) return false ; } return true ; } ll factor[100 ]; int tol;ll gcd (ll a, ll b) { if (a == 0 )return 1 ; if (a < 0 ) return gcd (-a, b); while (b) { ll t = a % b; a = b; b = t; } return a; } ll Pollard_rho (ll x, ll c) { ll i = 1 , k = 2 ; ll x0 = rand () % x; ll y = x0; while (1 ) { i++; x0 = (mult_mod (x0, x0, x) + c) % x; ll d = gcd (y - x0, x); if (d != 1 && d != x) return d; if (y == x0) return x; if (i == k) { y = x0; k += k; } } } void findfac (ll n) { if (Miller_Rabin (n)) { factor[tol++] = n; return ; } ll p = n; while (p >= n)p = Pollard_rho (p, rand () % (n - 1 ) + 1 ); findfac (p); findfac (n / p); } int main () { ll n; int t; cin >> t; while (t--) { cin >> n; if (n == 1 ) { continue ; } tol = 0 ; findfac (n); sort (factor, factor + tol); for (int i = 0 ; i < tol; i++) { cout << factor[i] << ' ' ; } } return 0 ; }
## 七, 乘法逆元
### 1.by EXGCD
O(logN)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 1 #include <bits/stdc++.h> 2 typedef long long ll; 3 ll exgcd (ll a,ll b,ll &x,ll &y) { 4 if (!b) { 5 x=1 ,y=0 ; 6 return a; 7 } 8 ll res=exgcd (b,a%b,y,x); 9 y-=a/b*x; 10 return res;11 }12 int main () 13 {14 ll a,p,x,y; 15 scanf ("%lld%lld" ,&a,&p);16 ll d=exgcd (a,p,x,y);17 printf (d==1 ?"%lld" :"-1" ,(x+p)%p);18 return 0 ;19 }
### 2.费马小定理
O(logP) 模为素数p的情况下
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 1 #include <bits/stdc++.h> 2 typedef long long ll; 3 ll quickpowmod (ll a,ll b,ll mod) { 4 ll ans=1 ; 5 while (b) { 6 if (b&1 ) ans=(ans*a)%mod; 7 b>>=1 ; 8 a=(a*a)%mod; 9 } 10 return ans;11 }12 int main () 13 {14 ll a,p; 15 scanf ("%lld%lld" ,&a,&p);16 ll inva=quickpowmod (a,p-2 ,p);17 printf ("%lld" ,inva);18 return 0 ;19 }
### 3.逆元打表
1 2 3 inv[1 ]=1 ; for (int i = 2 ; i <= n; i++) inv[i]=(long long )(p - p/i) * inv[p % i] % p;
## 八, 矩阵快速幂
模板是求解 : F(n)= F(n-1)+F(n-2)+F(n-3)
具体题目具体分析, 推算矩阵
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 #include <bits/stdc++.h> using namespace std;typedef long long ll;const ll MOD = 1000000007 ;struct mat { ll a[3 ][3 ]; }; mat mat_mul (mat a,mat b) { mat res; memset (res.a,0 ,sizeof (res.a)); for (int i = 0 ; i<=2 ; i++) { for (int j = 0 ; j<=2 ; j++) { for (int k = 0 ; k<=2 ; k++) { res.a[i][j] += ( a.a[i][k] * b.a[k][j])%MOD; } } } return res; } ll mat_pow (ll n) { mat c,res; memset (c.a,0 ,sizeof (c.a)); memset (res.a,0 ,sizeof (res.a)); c.a[0 ][0 ] = 1 ; c.a[0 ][1 ] = 1 ; c.a[0 ][2 ] = 1 ; c.a[1 ][0 ] = 1 ; c.a[2 ][1 ] = 1 ; for (int i=0 ; i<3 ; i++) { res.a[i][i] = 1 ; } while (n) { if (n&1 ) { res = mat_mul (res,c); } c = mat_mul (c,c); n >>= 1 ; } return (res.a[0 ][0 ])%MOD; } int main () { ll n; while (scanf ("%lld" ,&n)!=EOF) { ll ans = mat_pow (n); printf ("%lld\n" ,ans); } }
## 九, 卢卡斯定理
### 1.Lucas (模数p是质数)
给定$n,m,p (1 <= n,m,p <= 10^5)$ , 且p为素数
求 $C^n\_{n+m} mod P$
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 #include <iostream> using namespace std;typedef long long LL;const LL N=1e5 +2 ; LL a[N]; void init (LL p) { a[1 ]=1 ; for (int i=2 ;i<=p;++i)a[i]=a[i-1 ]*i%p; } void exgcd (LL a,LL b,LL &x,LL &y) { if (!b){ x=1 ; y=0 ; return ; } exgcd (b,a%b,y,x); y-=a/b*x; } LL ksm (LL x,LL n,LL mod) { LL ans=1 ; while (n){ if (n&1 )ans=ans*x%mod; n>>=1 ; x=x*x%mod; } return ans; } LL C (LL n,LL m,LL p) { if (n==m||m==0 )return 1 ; if (n<m)return 0 ; if (m*2 >n)m=n-m; return a[n]*ksm (a[m]*a[n-m],p-2 ,p)%p; } LL lucas (LL n,LL m,LL p) { if (!m)return 1 ; return lucas (n/p,m/p,p)*C (n%p,m%p,p)%p; } int main () { ios::sync_with_stdio (false ); LL T,n,m,p; cin>>T; while (T--){ cin>>n>>m>>p; init (p); cout<<lucas (n+m,m,p)<<endl; } return 0 ; }
### 2\. EX-Lucas (模数p不一定是质数的时候)
求$C^m\_n mod P$
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 #include <bits/stdc++.h> using namespace std;typedef long long LL;const LL N=1e5 +9 ;LL A[N],M[N]; LL ksm (LL x,LL n,LL mod) { LL ans=1 ; while (n){ if (n&1 )ans=ans*x%mod; n>>=1 ,x=x*x%mod; } return ans; } void exgcd (LL a,LL b,LL &x,LL &y) { if (!b)x=1 ,y=0 ; else exgcd (b,a%b,y,x),y-=a/b*x; } LL inv (LL a,LL p) { LL x,y; exgcd (a,p,x,y); return (x+p)%p?x:x+p; } LL get (LL n,LL pi,LL p) { if (!n)return 1 ; LL ans=1 ; if (n/p){ for (LL i=2 ;i<=p;++i)if (i%pi)ans=ans*i%p; ans=ksm (ans,n/p,p); } for (LL i=2 ;i<=n%p;++i)if (i%pi)ans=ans*i%p; return ans*get (n/pi,pi,p)%p; } LL exlucas (LL n,LL m,LL pi,LL p) { LL nn=get (n,pi,p); LL mm=get (m,pi,p); LL nm=get (n-m,pi,p); LL k=0 ; for (LL i=n;i;i/=pi)k+=i/pi; for (LL i=m;i;i/=pi)k-=i/pi; for (LL i=n-m;i;i/=pi)k-=i/pi; return nn*inv (mm,p)*inv (nm,p)*ksm (pi,k,p)%p; } LL crt (LL len,LL Lcm) { LL ans=0 ; for (LL i=1 ;i<=len;++i){ LL Mi=Lcm/M[i]; ans=((ans+A[i]*inv (Mi,M[i])*Mi)%Lcm+Lcm)%Lcm; } return ans; } int main () { ios::sync_with_stdio (false ); LL n,m,P,num; while (cin>>n>>m>>P){ if (n<m){ cout<<0 <<endl; continue ; } num=0 ; memset (A,0 ,sizeof (A)); memset (M,0 ,sizeof (M)); for (LL x=P,i=2 ;i<=P;++i) if (x%i==0 ){ M[++num]=1 ; while (x%i==0 ){ M[num]*=i; x/=i; } A[num]=exlucas (n,m,i,M[num])%P; } cout<<crt (num,P)<<endl; } return 0 ; }
## 十, 多项式算法
### FFT(根据复数, 答案是小数, 非精确)
### NTT
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 #include <bits/stdc++.h> #define LL long long int #define Redge(u) for (int k = h[u],to; k; k = ed[k].nxt) #define REP(i,n) for (int i = 1; i <= (n); i++) #define BUG(s,n) for (int i = 1; i <= (n); i++) cout<<s[i]<<' ' ; puts("" ); using namespace std;const int maxn = 4000005 , maxm = 100005 , INF = 1000000000 ;inline int read () { int out = 0 , flag = 1 ; char c = getchar (); while (c < 48 || c > 57 ) { if (c == '-' ) flag = -1 ; c = getchar (); } while (c >= 48 && c <= 57 ) { out = (out << 3 ) + (out << 1 ) + c - 48 ; c = getchar (); } return out * flag; } const int G = 3 , P = (119 << 23 ) + 1 ;int n, m, L, R[maxn];int A[maxn], B[maxn];int qpow (int a, int b) { int ans = 1 ; for (; b; b >>= 1 , a = 1ll * a * a % P) if (b & 1 ) ans = 1ll * ans * a % P; return ans; } void NTT (int * a, int f) { for (int i = 0 ; i < n; i++) if (i < R[i]) swap (a[i], a[R[i]]); for (int i = 1 ; i < n; i <<= 1 ) { int gn = qpow (G, (P - 1 ) / (i << 1 )); for (int j = 0 ; j < n; j += (i << 1 )) { int g = 1 ; for (int k = 0 ; k < i; k++, g = 1ll * g * gn % P) { int x = a[j + k], y = 1ll * g * a[j + k + i] % P; a[j + k] = (x + y) % P; a[j + k + i] = (x - y + P) % P; } } } if (f == 1 ) return ; int nv = qpow (n, P - 2 ); reverse (a + 1 , a + n); for (int i = 0 ; i < n; i++) a[i] = 1ll * a[i] * nv % P; } int main () { n = read (); m = read (); for (int i = 0 ; i <= n; i++) A[i] = read (); for (int i = 0 ; i <= m; i++) B[i] = read (); m = n + m; for (n = 1 ; n <= m; n <<= 1 ) L++; for (int i = 0 ; i < n; i++) R[i] = (R[i >> 1 ] >> 1 ) | ((i & 1 ) << (L - 1 )); NTT (A, 1 ); NTT (B, 1 ); for (int i = 0 ; i < n; i++) A[i] = 1ll * A[i] * B[i] % P; NTT (A, -1 ); for (int i = 0 ; i <= m; i++) printf ("%d " , A[i]); return 0 ; }
## 十一, 康托展开/逆康托展开
### 康托展开, 某排列映射一个值
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 const int fact[10 ] = {1 , 1 , 2 , 6 , 24 , 120 , 720 , 5040 , 40320 , 362880 };int contor (const vector<int >& permutation) { int num = 0 ; int len = permutation.size (); for (int i = 0 ; i < len; ++i) { int cnt = 0 ; for (int j = i + 1 ; j < len; ++j) if (permutation[i] > permutation[j]) ++cnt; num += cnt * fact[len - i - 1 ]; } return num + 1 ; }
### 反求排列
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 const int fact[10 ] = {1 , 1 , 2 , 6 , 24 , 120 , 720 , 5040 , 40320 , 362880 };vector<int > revContor (int bits, int num) { num = num - 1 ; vector<bool > vis (bits + 1 , false ) ; vector<int > permutation (bits, -1 ) ; int n, residue = num; for (int i = 0 ; i < bits; ++i) { n = residue / (fact[bits - i - 1 ]); residue = residue % (fact[bits - i - 1 ]); for (int j = 1 ; j <= bits; ++j) { if (!vis[j] && !(n--)) { vis[j] = true ; permutation[i] = j; break ; } } } return permutation; }
## 十二, 日期计算
### 吉姆拉尔森公式
1 2 3 4 5 6 7 8 int CaculateWeekDay (int y, int m, int d) { if (m == 1 || m == 2 ) { m += 12 ; y--; } int iWeek = (d + 2 * m + 3 * (m + 1 ) / 5 + y + y / 4 - y / 100 + y / 400 ) % 7 ; return iWeek; }
## 十三, 容斥原理
递归/暴力
若在n较大的情况下,应该采用递归进行计算,位运算计算循环次数太多,反而很慢。
在集合1,2,3…600,中,取出能被2,3,5,整除的数的个数
### 位操作:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 #include <bits/stdc++.h> using namespace std;int p[] = { 2 ,3 ,5 };int cnt = 3 ;int cal (int n = 600 ) { int res = 0 ; for (int i = 1 ; i < (1 << cnt); i++) { int t = i, tmp = 1 , k = 0 ; int len = 0 ; while (t) { if (t & 1 ) { tmp *= p[k]; len++; } t >>= 1 ; k++; } if (len & 1 ) res += n / tmp; else res -= n / tmp; } return res; } int main () { cout << cal () << endl; return 0 ; }
### 递归版
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 #include <bits/stdc++.h> using namespace std;int a[] = { 2 ,3 ,5 };int b = 600 ;int sum = 0 ;int n = 3 ;void dfs (int i, int num, int x, int mu) { if (num == x) { sum += b / mu; return ; } if (i == n) return ; dfs (i + 1 , num + 1 , x, mu * a[i]); dfs (i + 1 , num, x, mu); } int rong () { int s = 0 ; for (int i = 1 ; i <= n; i++) { sum = 0 ; dfs (0 , 0 , i, 1 ); if (i & 1 ) s += sum; else s -= sum; } return s; } void dfs (int i, int mu, int num) { mu *= a[i]; if (num % 2 ) sum += b / mu; else sum -= b / mu; for (int j = i + 1 ; j < 3 ; j++) { dfs (j, mu, num + 1 ); } } int main () { printf ("%d\n" , rong ()); sum = 0 ; for (int i = 0 ; i < 3 ; i++) { dfs (i, 1 , 1 ); } printf ("%d\n" , sum); return 0 ; }
### 从 m 中颜色中选择 k 种, 给 n 个盒子染色, 相邻盒子颜色不同
### 打表求 组合数(取模)
用不超过 k 种颜色涂完 n 个盒子的方案数为 $ x\_k=k(k−1)^{n−1} $
n 个盒子,用 m 种颜色中选出 正好 k 种颜色,且相邻盒子颜色不同的方案数为:
$$
ans = C\_{m}^{k} t m p=C\_{m}^{k} x\_{k}-C\_{m}^{k-1} x\_{k-1}+C\_{m}^{k-2} x\_{k-2}+\\ldots+(-1)^{k} x\_{1}
$$
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 #include<stdio.h> #include<cstdio> #include<iostream> using namespace std; #define ll long long const ll mod=1e9+7; const int maxn=1e6+10; int T; ll p[maxn], invp[maxn], inv[maxn]; //(a^b)%p ll qp(ll a, ll b){ a %= mod; ll ret = 1; while(b){ if (b & 1) ret = ret * a % mod; a = (a * a) % mod; b >>= 1; } return ret; } void init(){ p[0] = 1; for(int i = 1; i <= 1000000; i++) p[i] = p[i - 1] * i % mod; invp[1000000] = qp(p[1000000], mod - 2); for(int i = 1000000 - 1; i>=0; i--) invp[i] = invp[i + 1] * (i + 1) % mod; inv[1] = 1; for(int i = 2; i <= 1000000; i++) inv[i] = 1ll * (mod - mod / i) * inv[mod % i] % mod; } ll C(ll n, ll m){ return p[n] * invp[m] % mod * invp[n - m] % mod; } int main() { init(); scanf("%d",&T); while(T--) { ll n,m,k; scanf("%lld%lld%lld",&n,&m,&k); ll Cmk = 1; for(int i=1; i<=k; i++) Cmk = Cmk * (m - k + i) % mod * inv[i] % mod; ll tmp = 0; for(int i=0; i<k; i++) if (i % 2 == 0) tmp = (tmp + C(k, k-i) * (k - i) % mod * qp(k - i - 1, n - 1) % mod + mod) % mod; else tmp = (tmp - C(k, k-i) * (k - i) % mod * qp(k - i - 1, n - 1) % mod + mod) % mod; ll ans = Cmk * tmp % mod; printf("%lld\n", ans); } return 0; }
## 十四, 欧拉降幂
输入n, m
长度n序列, m是mod
q个查询, 每次输入l, r
求:
$$ w\[l\]^{w\[l+…\]^{w\[r\]}}$$
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 #include <bits/stdc++.h> using ll = long long ;using namespace std;#define MOD(a,b) a>=b?a%b+b:a const int maxn = 1e6 + 10 ;int is_prime[maxn];int prime[maxn];int phi[maxn];int cnt;void init () { cnt = 0 ; is_prime[1 ] = 1 ; phi[1 ] = 1 ; for (int i = 2 ; i <= maxn - 10 ; i++) { if (!is_prime[i]) { prime[++cnt] = i; phi[i] = i - 1 ; } for (int j = 1 ; j <= cnt; j++) { if (prime[j] * i >= maxn) break ; is_prime[i * prime[j]] = 1 ; if (i % prime[j] == 0 ) { phi[i * prime[j]] = phi[i] * prime[j]; break ; } else { phi[i * prime[j]] = phi[i] * (prime[j] - 1 ); } } } } unordered_map<int , int > ma; int n;ll w[maxn], m, l, r, x; ll qpow (ll a, ll b, ll mod) { ll res = 1 ; while (b) { if (b & 1 ) res = MOD (res * a, mod); a = MOD (a * a, mod); b >>= 1 ; } return res; } ll getphi (ll x) { if (x <= maxn - 10 ) { return phi[x]; } ll ans = x; if (ma.count (x)) return ma[x]; for (int i = 1 ; 1ll *prime[i] * prime[i] <= x; i++) { if (x % prime[i] == 0 ) { ans = ans / prime[i] * (prime[i] - 1 ); while (x % prime[i] == 0 ) { x /= prime[i]; } } } if (x > 1 ) ans = ans / x * (x - 1 ); ma[x] = ans; return ans; } ll solve (int ceng, ll mod) { if (ceng == r || mod == 1 ) return MOD (w[ceng], mod); return qpow (w[ceng], solve (ceng + 1 , getphi (mod)), mod); } int main () { init (); scanf ("%d %lld" , &n, &m); for (int i = 1 ; i <= n; i++) { scanf ("%d" , &w[i]); } int q; scanf ("%d" , &q); while (q--) { scanf ("%lld %lld" , &l, &r); cout << solve (l, m) % m << '\n' ; } return 0 ; }
## 十五, 高斯消元
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 const double EPS = 1E-9 ;int n;vector<vector<double > > a (n, vector <double >(n)); double det = 1 ;for (int i = 0 ; i < n; ++i) { int k = i; for (int j = i + 1 ; j < n; ++j) if (abs (a[j][i]) > abs (a[k][i])) k = j; if (abs (a[k][i]) < EPS) { det = 0 ; break ; } swap (a[i], a[k]); if (i != k) det = -det; det *= a[i][i]; for (int j = i + 1 ; j < n; ++j) a[i][j] /= a[i][i]; for (int j = 0 ; j < n; ++j) if (j != i && abs (a[j][i]) > EPS) for (int k = i + 1 ; k < n; ++k) a[j][k] -= a[i][k] * a[j][i]; } cout << det;
## 十六, 杜教BM
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 #include <bits/stdc++.h> using namespace std;typedef long long ll;typedef long long ll;typedef vector<int > VI;const int maxn = 10005 ;const ll mod = 1e9 + 7 ;const int INF = 0x3f3f3f3f ;const double eps = 1e-9 ;ll fast_mod (ll a, ll n, ll Mod) { ll ans = 1 ; a %= Mod; while (n) { if (n & 1 ) ans = (ans * a) % Mod; a = (a * a) % Mod; n>>= 1 ; } return ans; } namespace linear_seq{ ll res[maxn], base[maxn], num[maxn], md[maxn]; vector<int > vec; void mul (ll *a, ll *b, int k) { for (int i = 0 ; i <2 * k; i++) num[i] = 0 ; for (int i = 0 ; i < k; i++) { if (a[i]) for (int j = 0 ; j < k; j++) num[i + j] = (num[i + j] + a[i] * b[j]) % mod; } for (int i = 2 * k - 1 ; i>= k; i--) { if (num[i]) for (int j = 0 ; j <vec.size (); j++) num[i - k + vec[j]] = (num[i - k + vec[j]] - num[i] * md[vec[j]]) % mod; } for (int i = 0 ; i < k; i++) a[i] = num[i]; } ll solve (ll n, VI a, VI b) { ll ans = 0 , cnt = 0 ; int k = a.size (); assert (a.size () == b.size ()); for (int i = 0 ; i < k; i++) md[k - 1 - i] = -a[i]; md[k] = 1 ; vec.clear (); for (int i = 0 ; i < k; i++) if (md[i]) vec.push_back (i); for (int i = 0 ; i < k; i++) res[i] = base[i] = 0 ; res[0 ] = 1 ; while ((1LL << cnt) <= n) cnt++; for (int p = cnt; p>= 0 ; p--) { mul (res, res, k); if ((n>> p) & 1 ) { for (int i = k - 1 ; i>= 0 ; i--) res[i + 1 ] = res[i]; res[0 ] = 0 ; for (int j = 0 ; j < vec.size (); j++) res[vec[j]] = (res[vec[j]] - res[k] * md[vec[j]]) % mod; } } for (int i = 0 ; i < k; i++) ans = (ans + res[i] * b[i]) % mod; if (ans < 0 ) ans += mod; return ans; } VI BM (VI s) { VI B (1 , 1 ) , C (1 , 1 ) ; int L = 0 , m = 1 , b = 1 ; for (int i = 0 ; i < s.size (); i++) { ll d = 0 ; for (int j = 0 ; j < L + 1 ; j++) d = (d + (ll)C[j] * s[i - j]) % mod; if (d == 0 ) m++; else if (2 * L <= i) { VI T = C; ll c = mod - d * fast_mod (b, mod - 2 , mod) % mod; while (C.size () < B.size () + m) C.push_back (0 ); for (int j = 0 ; j < B.size (); j++) C[j + m] = (C[j + m] + c * B[j]) % mod; L = i + 1 - L, B = T, b = d, m = 1 ; } else { ll c = mod - d * fast_mod (b, mod - 2 , mod) % mod; while (C.size () < B.size () + m) C.push_back (0 ); for (int j = 0 ; j < B.size (); j++) C[j + m] = (C[j + m] + c * B[j]) % mod; m++; } } return C; } int gao (VI a, ll n) { VI c = BM (a); c.erase (c.begin ()); for (int i = 0 ; i < c.size (); i++) c[i] = (mod - c[i]) % mod; return solve (n, c, VI (a.begin (), a.begin () + c.size ())); } } int main () { ll t, n; scanf ("%d" , &t); while (t--) { scanf ("%lld" , &n); printf ("%lld\n" , linear_seq::gao (VI{3 , 9 , 20 , 46 , 106 , 244 , 560 , 1286 , 2956 , 6794 }, n-1 )); } return 0 ; }
## 十七, 拉格朗日插值
高数的应用: 根据导数, 控制函数大小
这里的应用: 给出n个点P(xi, yi),将过这个点的最多n - 1次的多项式记为f(x),求f(k)的值
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 #include <cstdio> const int maxn = 2010 ;using ll = long long ;ll mod = 998244353 ; ll n, k, x[maxn], y[maxn], ans, s1, s2; ll powmod (ll x, ll n) { ll ret = 1ll ; while (n) { if (n & 1 ) ret = ret * x % mod; x = x * x % mod; n >>= 1 ; } return ret; } ll inv (ll x) { return powmod (x, mod - 2 ); }int main () { scanf ("%lld%lld" , &n, &k); for (int i = 1 ; i <= n; i++) scanf ("%lld%lld" , x[i], y[i]); for (int i = 1 ; i <= n; i++) { s1 = y[i] % mod; s2 = 1ll ; for (int j = 1 ; j <= n; j++) if (i != j) s1 = s1 * (k - x[j]) % mod, s2 = s2 * (x[i] - x[j]) % mod; ans += s1 * inv (s2) % mod; } printf ("%lld\n" , (ans % mod + mod) % mod); return 0 ; }
## 十八, Catalan 数列
常常用来求解这样一类问题:
1. 有2n个人排成一行进入剧场。入场费5元。其中只有n个人有一张5元钞票,另外n人只有10元钞票,剧院无其它钞票,问有多少中方法使得只要有 10 元的人买票,售票处就有 5 元的钞票找零?
2. 在圆上选择2n个点,将这些点成对连接起来使得所得到的n条线段不相交的方法数?
3. n个不同的数依次进栈,求不同的出栈结果的种数?
这一类, 有两种选择, 每种选择只能选择有限次数, 并且两种选择的可选择次数相同, 问多少种方案的问题.
1 2 3 4 5 6 7 8 9 10 11 12 #include <iostream> using namespace std;int n;long long f[25 ];int main () { f[0 ] = 1 ; cin >> n; for (int i = 1 ; i <= n; i++) f[i] = f[i - 1 ] * (4 * i - 2 ) / (i + 1 ); cout << f[n] << endl; return 0 ; }
## 十九, 常见的期望和方差
#### 0-1分布
1 次实验
P{x = 0} = p
P{x = 1} = 1 - p
E(x) = p;
D(x) = pq;
#### 二项分布 B(n, p)
n次独立重复实验
P{x = k} = p
E(x) = np;
D(x) = npq;
#### 几何分布
在n次独立重复实验中, 前k-1次皆失败,第k次成功的概率
P{x = k} = p
E(x) = 1/p;
D(x) = (1 - p) / (p \* p);