zhugezy的完美数学教室(v2.0)
你找到了这个菜鸡的板子。。。
2019/09/20:马上要ccpc了,更新一下,过十几个小时应该还会加一点东西(?)
zhugezy的完美数学教室
一些结论
\(2C_n^2 + n = n^2,C_n^m=C_{n - 1}^m+C_{n-1}^{m-1}\)
\(\sum_{i = 1}^niC_n^i=n2^{n-1},\sum_{i=1}^ni^2C_n^i=n(n+1)2^{n-2}\)
\(C_n^{x+1}+2\sum_{i=0}^xC_n^i=\sum_{i=0}^{x+1}C_{n+1}^i\)
\(\sum_{i=1}^n(-1)^{i-1}\frac{C_n^i}{i}=\sum_{i=1}^n\frac{1}{i}\)
\(\sum_{i=0}^n(C_n^i)^2=C_{2n}^n\)
1~n与n互质的数的和:\(\frac{n*\phi(n)+[n==1]}{2}\)
\(\sum_{i=1}^ngcd(i,n)=\sum_{d|n}d\phi(n/d)\)
1p模p的所有逆元值对应1p中所有的数,例如p=7时,1~6分别对应1,4,5,2,3,6
$1 \leq i \leq n, \lfloor\frac{n}{i}\rfloor $最多只有\(2\sqrt{n}\)种取值
\(\lfloor {\frac{n}{\lfloor{\frac{n}{i}}\rfloor}}\rfloor\)是\(\lfloor\frac{n}{i}\rfloor\)取值相同的一段区间的右端点
\(\lfloor\frac{n}{ab}\rfloor=\lfloor\frac{\lfloor\frac{n}{a}\rfloor}{b}\rfloor=\lfloor\frac{\lfloor\frac{n}{b}\rfloor}{a}\rfloor\)(对\(\lceil\rceil\)也成立)
n的所有因子之和(对n的质因数分解式)\(\sigma(n) = \frac{p^{e_1 + 1}_1 - 1}{p_1-1} *\frac{p^{e_2 + 1}_2 - 1}{p_2-1}*...*\frac{p^{e_k + 1}_k - 1}{p_k-1}\)
n个数的全排列的逆序对数和:\(n!\frac{n(n-1)}{4}\)
\(Cayley\)公式:\(n\)个带编号的结点组成\(s\)棵树的森林,其中指定\(s\)个不同的结点属于不同的树,这样的组成方案一共有\(F(n,s)=sn^{n-s-1}\)种。
特别地,\(n\)个结点带编号的完全图一共有\(n^{n-2}\)个生成树,即\(n\)个结点带编号的无根树有\(n^{n-2}\)个。
威尔逊定理:\(p\)为质数\(\ \ iff \ \ p|(p-1)!+1\)
\(gcd(a^n−1,a^m−1)=a^{gcd(n,m)}−1\)
\(a>b,gcd(a,b)=1\)时有\(gcd(a^m-b^m,a^n-b^n)=a^{gcd(m,n)}-b^{gcd(m,n)}\)
设\(G(n)=gcd(C_n^1,C_n^2,...,C_n^{n-1}),\)则对任意质数\(p\)和自然数\(i\),有\(G(p^i)=p\),否则\(G(n)=1\).
\((n+1)lcm(C_n^0,C_n^1,...,C_n^n)=lcm(1,2,...,n+1)\)
对质数\(p,(x_1+x_2+...+x_n)^p\equiv(x_1^p+x_2^p+...+x_n^p)\ (mod\ p)\).
斐波那契数列
\(F_i=\frac{1}{\sqrt5}[(\frac{1+\sqrt5}{2})^i-(\frac{1-\sqrt5}{2})^i]\)
\(gcd(F_i,F_j) = F_{gcd(i,j)}\)
\(F_n^2-F_{n-1}F_{n+1}=(-1)^{n-1}\)
\(F_1+F_3+...+F_{2n-1}=F_{2n}-F_2+F_1\)
\(F_2+F_4+...+F_{2n}=F_{2n+1}-F_1\)
\(F_1^2+F_2^2+...+F_n^2=F_nF_{n+1}\)
\(F_{2n-2m-2}(F_{2n}+F_{2n+2})=F_{2m+2}+F_{4n-2m}(n > m \geq-1,n\geq1)\)
\(\frac{F_{2n}}{F_n}=F_{n-1}+F_{n+1}\)
中国剩余定理
且\(m_1, m_2, ...,m_n\)两两互质,则对任意的\(a_1, a_2, ...,a_n\)同余方程组有解,
设\(M = m_1*m_2*...*m_n\) 是整数\(m_1, m_2, ...,m_n\)的乘积,并设 \(M_i\)是除了\(m_i\)以外的n-1个整数的乘积。
设 \(t_i\)为\(M_i\) 模 \(m_i\) 意义下的逆元)
\(M_i * t_i\equiv{1}\pmod{m_i}(1\leq i \leq n)\)
方程组 (S) 的通解形式为
$x = kM + ∑(a_it_i*M_i),i∈{1,2,3,…,n},k∈Z $
在模 M 的意义下,方程组(S)只有一个解如下图所示
//求M%A=a,M%B=b,...中的M,其中A,B,C...互质
int CRT(int a[],int m[],int n){
int M = 1;
int ans = 0;
for(int i=1; i<=n; i++)
M *= m[i];
for(int i=1; i<=n; i++){
int x, y;
int Mi = M / m[i];
ex_gcd(Mi, m[i], x, y);
ans = (ans + Mi * x * a[i]) % M;
}
if(ans < 0) ans += M;
return ans;
}
#include<bits/stdc++.h>/*模不互质*/
using namespace std;
int a[10],m[10],tong[103];
int exgcd(int a,int &x,int b,int &y)
{
if(b==0){x=1;y=0;return a;}
int x2,y2;
int gcd=exgcd(b,x2,a%b,y2);
x=y2;
y=x2-a/b*y2;
return gcd;
}
int main()
{
freopen("crt.in","r",stdin);
freopen("crt.out","w",stdout);
int T;
scanf("%d",&T);
while(T--)
{
int n;
scanf("%d",&n);
scanf("%d%d",&a[1],&m[1]);
//用lcm求当前到达过的方程中,所有模数的lcm,sum是最后的结果,但是每一次都会用
//解出的x去更新它,fail判断是否有解
int lcm=m[1],sum=a[1],fail=0;
for(int i=2;i<=n;i++)
{
int x,y;
scanf("%d%d",&a[i],&m[i]);
a[i]=((a[i]-sum)%m[i]+m[i])%m[i];//sum就是要求的x,不断更新
int d=exgcd(lcm,x,m[i],y);//lcm*x+m[i]*y+sum=a[i]才能保证x%m[i]=a[i]
if(a[i]%d==0) x=x*(a[i]/d)%m[i];
else fail=1;
sum+=x*lcm;//注意是乘lcm哦!不然可能会出现前面的方程冲突
lcm=lcm/d*m[i];//到现在这一个方程了,所有模数的最小公倍数
sum=(sum%lcm+lcm)%lcm;
}
if(fail) printf("No\n");
else printf("%d\n",sum);
}
}
lucas定理 求 \(C_n^k\% p\)(p为素数,p规模小于1e5)
表达式:C(n,m)%p=C(n/p,m/p)*C(n%p,m%p)%p。(可以递归)
递归方程:(C(n%p, m%p)*Lucas(n/p, m/p))%p。(递归出口为m==0,return 1)
扩展lucas
输入三个整数\(n,k,p\),输出\(C_n^k\%p\).\(1\leq k \leq n \leq 1e18, 2 \leq p \leq 1e6\),\(p\)不一定是质数。
#include<bits/stdc++.h>
using namespace std;
const int maxn = 100000 + 10;
typedef long long LL;
LL Pow(LL n, LL m, LL mod) {
LL ans = 1;
while(m > 0) {
if(m & 1) ans = (LL)ans * n % mod;
n = (LL)n * n % mod; m >>= 1;
}
return ans;
}
LL Pow(LL n,LL m) {
LL ans = 1;
while(m > 0) {
if(m & 1) ans = ans * n;
n = n * n; m >>= 1;
}
return ans;
}
LL x, y;
LL exgcd(LL a, LL b) {
if(a == 0) {
x = 0, y = 1;
return b;
}LL r = exgcd(b%a, a);
LL t = x; x = y - (b/a)*x; y = t;
return r;
}
LL rev(LL a, LL b) { exgcd(a, b); return ((x % b) + b) % b; }
LL Calc(LL n, LL p, LL t) {
if(n == 0) return 1;
LL s = Pow(p, t), k = n / s, tmp = 1;
for(LL i=1; i<=s; i++) if(i % p) tmp = (LL)tmp * i % s;
LL ans = Pow(tmp, k, s);
for(LL i=s*k + 1; i<=n; i++) if(i % p) ans = (LL)ans * i % s;
return (LL)ans * Calc(n / p, p, t) % s;
}
LL C(LL n, LL m, LL p, LL t) {
LL s = Pow(p, t), q = 0;
for(LL i=n; i; i/=p) q += i / p;
for(LL i=m; i; i/=p) q -= i / p;
for(LL i=n-m; i; i/=p) q -= i / p;
LL ans = Pow(p, q);
LL a = Calc(n, p, t), b = Calc(m, p, t), c = Calc(n-m, p, t);
return (LL)(ans * a * rev(b, s) * rev(c, s)) % s;
}
LL China(LL A[], LL M[], LL cnt) {
LL ans = 0, m, n = 1;
for(LL i=1; i<=cnt; i++) n *= M[i];
for(LL i=1; i<=cnt; i++) {
m = n / M[i];
exgcd(M[i], m);
ans = (ans + (LL)y * m * A[i]) % n;
}
return (ans + n) % n;
}
LL A[maxn], M[maxn], cnt;
LL Lucas(LL n, LL m, LL mod) {
for(LL i=2; i*i <= mod; i++) if(mod % i == 0) {
LL t = 0;
while(mod % i == 0) t++, mod /= i;
M[++cnt] = Pow(i, t);
A[cnt] = C(n, m, i, t);
}if(mod > 1) {
M[++cnt] = mod;
A[cnt] = C(n, m, mod, 1);
}
return China(A, M, cnt);
}
LL n, k, p;
int main() {
#ifndef ONLINE_JUDGE
freopen("input.in", "r", stdin);
freopen("output.out", "w", stdout);
#endif
cin >> n >> k >> p;
cout << Lucas(n, k, p) << endl;
return 0;
}
一些常见的组合数学数
卡特兰数\(C_n\)
1, 2, 5, 14, 42,132, 429, 1430, 4862, 16796,58786, 208012, 742900, 2674440, 9694845, 35357670, 129644790, 477638700, 1767263190, 6564120420, 24466267020, 91482563640, 343059613650, 1289904147324, 4861946401452
通项\(C_n = \frac{1}{n + 1} C^{n}_{2n}\)
性质\(C_n = C^{n}_{2n} - C^{n - 1}_{2n}\) \(C_0 = 1, C_{n+1} = \frac{2(2n+1)}{n+2}C_n\)
\(C_{n+1} = \sum_{i = 0}^nC_iC_{n - i}\) 若\(C_n\)是奇数,那么\(n=2^k-1\)
斯特林数
第一类Stirling数
\(s(n,m)\)表示将 n 个不同元素构成m个圆排列的数目。
无符号第一类斯特林数\(S_u(n,m)\) 带符号第一类斯特林数\(S_s(n, m)\)
\(x^{n\downarrow}=\prod_{i=0}^{n-1}(x-i)=\sum_{k=0}^ns_s(n,k)x^k\)
\(x^{n\uparrow}=\prod_{i=0}^{n-1}(x+i)=\sum_{k=0}^ns_u(n,k)x^k\)
有\(s_s(n,m)=(-1)^{n+m}s_u(n,m)\).
有\(s_u(n+1,m)=s_u(n,m-1)+ns_u(n,m)\)
\(s_s(n+1,m)=s_s(n,m-1)-ns_s(n,m)\)
\(s_u(0,0)=1,s_u(n,0)=0(n \geq 1), s_u(n,n)=1,s_u(n,1)=(n-1)!\)
\(s_u(n,n-1)=C_n^2,s_u(n,2)=(n-1)!*\sum_{i=1}^{n-1}\frac{1}{i},s_u(n,n-2)=2C_n^3+3C_n^4\)
\(\sum_{k=0}^ns_s(n,k)=0^n\),注意\(0^0=1\).
第二类Stirling数
\(n\)个不同元素拆分成\(m\)个无区别的集合的方案数
\(S(n,m)=\frac{1}{m!}\sum_{k=0}^m(-1)^kC_m^k(m-k)^n\)
\(S(n+1,m)=S(n,m-1)+mS(n,m)\)
两类数的转化\(\sum_{k=0}^nS(n,k)s(k,m)=\sum_{k=0}^ns(n,k)S(k,m)\)
贝尔数
集合划分的数目,如\(B_3=5,\)表示\(\{a,b,c\}\)的五个划分为\(\{\{a\},\{b\}, \{c\}\},\{\{a\}, \{b, c\}\},\{\{b\}, \{a, c\}\},\{\{c\}, \{a, b\}\},\{\{a,b,c\}\}\)
\(B_0=B_1=1\)
$1, 1, 2, 5, 15, 52, 203, 877, 4140, 21147, 115975, 678570, 4213597, 27644437, 190899322, 1382958545, 10480142147, 82864869804, 682076806159, 5832742205057, $$51724158235372, 474869816156751, 4506715738447323, 44152005855084346, 445958869294805289, 4638590332229999353, 49631246523618756274......$
\(B_{n+1}=\sum\limits_{k=0}^nC_n^kB_k\)
\(B_n=\frac{1}{e}\sum\limits_{k=0}^{\infin}\frac{k^n}{k!}\)
\(Touchard\)同余:若\(p\)为质数,则\(B_{p+n}\equiv B_n + B_{n+1}(mod\ p)\).
\(B_n=\sum\limits_{k=1}^nS(n,k)\),\(S\)是第二类Stirling数.
常用数论函数
欧拉函数
\(\phi(p) = p - 1, \phi(p^k) = (p - 1)*p^{k - 1}\)
\(时,gcd(m,n) = 1时, \phi(mn) = \phi(m)\phi(n)\)
\(\sum_{d|n}\phi(d) = n\)
int getphi(int n)
{
int ans = n;
for (int i = 2; i * i <= n; ++i)
{
if (n % i == 0)
{
ans -= ans / i;
while (n % i == 0)
n /= i;
}
}
if (n > 1)
ans -= ans / n;
return ans;
}
欧拉定理(n为素数即为费马小定理)
\(则有gcd(a,n) = 1,则有a^{\phi(n)}\equiv{1}\mod{n}\)
推广:欧拉降幂
LL Mod(LL b, LL mod){return b < mod?b:b%mod+mod;}
莫比乌斯函数
n质因数分解后,\(如果有幂次大于的为质因子个数\mu(1) = 1;\mu(n) = 0如果有幂次大于1的;\mu(n) = (-1)^k,k为质因子个数\)
n>2时,n的所有约数对应的莫比乌斯函数值之和为0 即\(\sum_{d|n}\mu(d) = 0\)
\(\sum_{d|n}\frac{\mu(d)}{d} = \frac{\phi(n)}{n}\)
等价于\(\phi(n)=\sum_{d|n}\frac{n}{d}\mu(d)=\sum_{d|n}d\mu(\frac{n}{d})\)
\([gcd(a,b) == 1] = \sum_{d|gcd(a,b)}\mu(d)=\sum_{d|a,d|b}\mu(d)\) 容斥
例题:求1-n中有多少个数和x互质
\(\sum_{d|x}\mu(d)*\lfloor\frac{n}{d}\rfloor\)容斥系数即为莫比乌斯函数
int mobius(int x)
{
int res=1;
int i,j;
for(i=2;i*i<=x;i++)
if(x%i==0)
{
x/=i;
if(x%i==0) mu[i]=0;
else mu[i]*=-1;
while(x%i==0) x/=i;
}
return res;
}
exgcd
#include <iostream>
#define LL long long
using namespace std;
LL exgcd(LL a, LL b, LL &p, LL &q)
{
if (a == 0 && b == 0) return -1;
if (b == 0)
{
p = 1; q = 0; return a;
}
LL d = exgcd(b, a % b, q, p);
q -= a / b * p;
return d;
}
逆元
1.费马小定理----mod为素数
qp(a, mod - 2);
2.exgcd
给定模数m,求a的逆相当于求解ax=1(mod m),这个方程可以转化为ax-my=1 , 然后套用求二元一次方程的方法,用扩展欧几里得算法求得一组x0,y0和gcd
检查gcd是否为1 ,gcd不为1则说明逆元不存在
若为1,则调整x0到0~m-1的范围中即可
#define LL long long
void extgcd(LL a,LL b,LL &d,LL &x,LL &y)
{
if(!b){ d=a; x=1; y=0;}
else{ extgcd(b,a%b,d,y,x); y-=x*(a/b); }
}
LL inv(LL a,LL n)
{
LL d,x,y;
extgcd(a,n,d,x,y);
return d==1?(x+n)%n:-1;
}
3.b|a
4.欧拉定理(m可以不为质数,但要a和m互质)
5.逆元打表
#define LL long long
const int N = 1e5 + 5;
int inv[N];
void getinv(int n, int p)
{
inv[1] = 1;
for (int i = 2; i <= n; ++i)
{
inv[i] = (LL) (p - p / i) * inv[p%i] % p;
}
}
线性筛(素数筛,欧拉函数,莫比乌斯函数,约数个数,最小因子)
#define MAXN 1000000
int prime[MAXN + 8];
bool notprime[MAXN + 8];
//int phi[MAXN + 8];
///int mu[MAXN + 8];
////int facnum[MAXN + 8], d[MAXN + 8];
int ind = 0;
void getprime()
{
//phi[1] = 1;
///mu[1] = 1;
////facnum[1] = 1;
for (int i = 2; i <= MAXN; ++i)
{
if(!notprime[i])
{
prime[++ind] = i;
//phi[i] = i - 1;
///mu[i] = -1;
////facnum[i] = 2, d[i] = 1;
/////minfactor[i] = i;
}
for (int j = 1; j <= ind && i * prime[j] <= MAXN; ++j)
{
notprime[i * prime[j]] = 1;
if(i % prime[j] == 0)
{
//phi[i * prime[j]] = phi[i] * prime[j];
///mu[i * prime[j]] = 0;
////facnum[i * prime[j]] = facnum[i] / (d[i] + 1) * (d[i] + 2), d[i * prime[j]] = d[i] + 1;
break;
}
//phi[i * prime[j]] = phi[i] * (prime[j] - 1);
///mu[i * prime[j]] = -mu[i];
////facnum[i * prime[j]] = facnum[i] * 2, d[i * prime[j]] = 1;
/////minfactor[i*prime[j]] = prime[j];
}
}
}
例子:前n个数的约数个数筛
facnum[i]表示i的约数个数
通过素数筛得到前n个数的约数个数非常巧妙,首先根据约数个数定理:
对于一个大于1正整数n可以分解质因数:\(n = p_1^{d_1} * p_2^{d_2} *... * p_k^{d_k}\),其中pi为素数
则n的正约数的个数就是
:facnum[n] = (1 + d1) * (1 + d2) * ... * (1 + dk)
我们需要一个辅助数组d[i],表示i的最小质因子的次幂,(最小的原因是素数筛里每次都是用最小的质因子来筛合数的),还是三种情况:
1.当i为素数时,facnum[i] = 2;d[i] = 1,很好理解
2.当i % p[j] != 0时,gcd(i, p[j]) =1,由积性函数的性质可得facnum[i * p[j]] = facnum[i] * facnum[p[j]] = facnum[i] * 2
d[i * p[j]] = 1(无平方因子)
3.当i % p[j] == 0时,出现平方因子,最小质因子的次幂加1,因此有
facnum[i * p[j]] = facnum[i] / (d[i] + 1) * (d[i] + 2)
d[i * p[j]] = d[i] + 1
康托展开
康托展开
康托展开是一个全排列到一个自然数的双射,常用于构建hash表时的空间压缩。设有n个数(1,2,3,4,…,n),可以有组成不同(n!种)的排列组合,康托展开表示的就是是当前排列组合在n个不同元素的全排列中的名次。
在(1,2,3,4,5)5个数的排列组合中,计算 34152的康托展开值。
首位是3,则小于3的数有两个,为1和2,a[5]=2,则首位小于3的所有排列组合为 a[0]*(5-1)!
第二位是4,则小于4的数有两个,为1和2,注意这里3并不能算,因为3已经在第一位,所以其实计算的是在第二位之后小于4的个数。因此a[4]=2
第三位是1,则在其之后小于1的数有0个,所以a[3]=0
第四位是5,则在其之后小于5的数有1个,为2,所以a[2]=1
最后一位就不用计算啦,因为在它之后已经没有数了,所以a[1]固定为0
根据公式:
X = 2 * 4! + 2 * 3! + 0 * 2! + 1 * 1! + 0 * 0! = 2 * 24 + 2 * 6 + 1 = 61
所以比 34152 小的组合有61个,即34152是排第62。
static const int FAC[] = {1, 1, 2, 6, 24, 120, 720, 5040, 40320, 362880}; // 阶乘
int cantor(int *a, int n)
{
int x = 0;
for (int i = 0; i < n; ++i) {
int smaller = 0; // 在当前位之后小于其的个数
for (int j = i + 1; j < n; ++j) {
if (a[j] < a[i])
smaller++;
}
x += FAC[n - i - 1] * smaller; // 康托展开累加
}
return x; // 康托展开值
}
逆康托展开
一开始已经提过了,康托展开是一个全排列到一个自然数的双射,因此是可逆的。即对于上述例子,在(1,2,3,4,5)给出61可以算出起排列组合为 34152。由上述的计算过程可以容易的逆推回来,具体过程如下:
用 61 / 4! = 2余13,说明a[5]=2,说明比首位小的数有2个,所以首位为3。
用 13 / 3! = 2余1,说明a[4]=2,说明在第二位之后小于第二位的数有2个,所以第二位为4。
用 1 / 2! = 0余1,说明a[3]=0,说明在第三位之后没有小于第三位的数,所以第三位为1。
用 1 / 1! = 1余0,说明a[2]=1,说明在第二位之后小于第四位的数有1个,所以第四位为5。
最后一位自然就是剩下的数2啦。
通过以上分析,所求排列组合为 34152。
static const int FAC[] = {1, 1, 2, 6, 24, 120, 720, 5040, 40320, 362880}; // 阶乘
void decantor(int x, int n)
{
vector<int> v; // 存放当前可选数
vector<int> a; // 所求排列组合
for(int i=1;i<=n;i++)
v.push_back(i);
for(int i=m;i>=1;i--)
{
int r = x % FAC[i-1];
int t = x / FAC[i-1];
x = r;
sort(v.begin(),v.end());// 从小到大排序
a.push_back(v[t]); // 剩余数里第t+1个数为当前位
v.erase(v.begin()+t); // 移除选做当前位的数
}
}
莫比乌斯反演+数论分块
如果存在函数F(x)和f(x),满足 \(F(n) = \sum_{d|n}f(d)\),那么
或如果\(F(n) = \sum_{n|d}f(d)\),那么
例题:求对于区间[a,b]内的整数x和[c, d]内的y,满足gcd(x,y)=k的数对的个数。
\(的对数F(t):gcd(x, y) \% t == 0 的(x,y)对数\) \(的对数f(t):gcd(x, y) == t 的(x, y)对数\)
考虑 $x \in [1,b] $ and $ y \in [1,d]$ : \(F(t) = \lfloor\frac{b}{t}\rfloor\lfloor\frac{d}{t}\rfloor\)
$Ans = f(k), gcd(x,y) = k \leftrightarrow gcd(x/k, y/k) = 1 \leftrightarrow $
\(gcd(x,y) = 1(x\in[1,b/k]\) and \(y\in[1,d/k])\)
因此直接把b,d都除上k,\(f(1) = \sum^{min(b/k,d/k)}_{i = 1}\mu(i)F(i)\)就是结果
当t接近最大值时,对一段连续的t它们的\(F(t)\)值都相同,因此可以打出莫比乌斯前缀和然后分块求解。
#include <bits/stdc++.h>
#define LL long long
#define pb push_back
#define MAXN 50000
using namespace std;
bool check[MAXN + 8];
int prime[MAXN + 8];
int mu[MAXN + 8];
int sum[MAXN + 8];
void Mobius()
{
memset(check, false, sizeof(check));
mu[1] = 1;
int tot = 0;
for (int i = 2; i <= MAXN; ++i)
{
if (!check[i])
{
prime[++tot] = i;
mu[i] = -1;
}
for (int j = 1; j <= tot; ++j)
{
if (i * prime[j] > MAXN)
break;
check[i * prime[j]] = true;
if (i % prime[j] == 0)
{
mu[i * prime[j]] = 0;
break;
}
else
{
mu[i * prime[j]] = -mu[i];
}
}
}
for (int i = 1; i <= MAXN; ++i)
{
sum[i] = sum[i - 1] + mu[i];
}
}
LL calc(int a, int b)
{
LL ret = 0;
int last = 0;
for (int i = 1; i <= min(a, b); i = last + 1)
{
last = min(a / (a / i), b / (b / i));
ret += 1LL * (sum[last] - sum[i - 1]) * (a / i) * (b / i);
}
return ret;
}
int main()
{
Mobius();
int T, a, b, c, d, k;
cin >> T;
for (int t = 1; t <= T; ++t)
{
cin >> a >> b >> c >> d >> k;
a--, c--;
cout << calc(b / k, d / k) - calc(a / k, d / k) - calc(b / k, c / k) + calc(a / k, c / k);//容斥
if (t != T) cout << endl;
}
return 0;
}
卷积
\(id(n)=n,\epsilon(n)=[n=1],I(n)=1,d(n)=\sum_{d|n}1,\sigma(n)=\sum_{d|n}d,\lambda(n)=(-1)^k\)
\(\mu*I=\epsilon,\phi*I=id,\phi=id*\mu,d=I*I,1=\mu*d,\sigma*\mu=id,\sigma=id*I\)
\(\epsilon\)卷任何\(f\)都为\(f\)
常见求和公式及其推导
1.\(\sum_{a=1}^n\sum_{b=1}^mgcd(a,b)=\sum_{d=1}^{min(n,m)}\phi(d)\lfloor\frac{n}{d}\rfloor\lfloor\frac{m}{d}\rfloor\)
不妨设\(n\leq m\),枚举\(d=gcd(a,b)\),反演得到左边=\(\sum_{d=1}^n\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}d\mu(i)\lfloor\frac{n}{id}\rfloor\lfloor\frac{m}{id}\rfloor\),变形,用\(t\)替换\(i*d\),得\(\sum_{t=1}^n\lfloor\frac{n}{t}\rfloor\lfloor\frac{m}{t}\rfloor\sum_{d|t}d\mu(\frac{t}{d})\),后面那个式子是典型卷积,所以原公式成立。
2.\(\sum_{a=1}^n\sum_{b=1}^mlcm(a,b)=\frac{1}{4}\sum_{d=1}^{min(n,m)}d\lfloor\frac{n}{d}\rfloor\lfloor\frac{m}{d}\rfloor(\lfloor\frac{n}{d}\rfloor+1)(\lfloor\frac{m}{d}\rfloor+1)\sum_{i|d}i\mu(i)\)
反演:设\(F(t)=\sum_{i=1}^n\sum_{j=1}^mij[t|gcd(i,j)],f(t)=\sum_{i=1}^n\sum_{j=1}^mij[t=gcd(i,j)]\),则有\(F(t)=\sum_{i=1}^n\sum_{j=1}^mij[t|i][t|j]=\frac{\lfloor\frac{n}{t}\rfloor\lfloor\frac{m}{t}\rfloor(\lfloor\frac{n}{t}\rfloor+1)(\lfloor\frac{m}{t}\rfloor+1)}{4}\).因此答案=\(\sum_{d=1}^nd\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\mu(i)i^2\frac{\lfloor\frac{n}{id}\rfloor\lfloor\frac{m}{id}\rfloor(\lfloor\frac{n}{id}\rfloor+1)(\lfloor\frac{m}{id}\rfloor+1)}{4}\).\(t=id\)来替换\(d\),套路求解。记\(G(d)=\sum_{i|d}i\mu(i)\),显然是积性函数,对质数\(p\)有\(G(p^k)=1-p\),随便推一推就能用线性筛预处理,前面再一分块,就能\(O(\sqrt n)\)求解。
3.\(\sum_{a=1}^n\sum_{b=1}^nlcm(a,b)=\sum_{i=1}^n(-i+2\sum_{j=1}^ilcm(i,j))\).
预处理,O(1)输出。
4.\(\sum_{i=1}^ngcd(i,n)=\sum_{d|n}\frac{n}{d}\phi(d)\).
枚举\(d\),老套路了。
5.\(\sum_{i=1}^n\sigma(i)=\sum_{i=1}^ni\lfloor\frac{n}{i}\rfloor\).
根据\(\sigma\)的定义轻松得出。左边=\(\sum_{i=1}^n\sum_{j|i}j=\sum_{j=1}^nj\sum_{i=1}^{\lfloor\frac{n}{j}\rfloor}1\).
一类数论函数的更快求和方法
杜教筛---\(O(n^{2/3})\)计算一类积性函数的前缀和
我们要求\(\sum_{i=1}^nf(i)\),其中\(f\)为积性函数。
构造两个积性函数\(g,h\)满足\(h=f*g\)
\(\sum_{i=1}^nh(i)=\sum_{i=1}^n\sum_{d|i}g(d)f(i/d)=\sum_{d=1}^ng(d)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}f(i)=\sum_{d=1}^ng(d)S_f(\lfloor\frac{n}{d}\rfloor)\)
因此\(g(1)S_f(n)=\sum_{i=1}^{n}h(i)-\sum_{i=2}^ng(d)S_f(\lfloor\frac{n}{d}\rfloor)\)
包含\(\phi(i),\mu(i),i\phi(i),i\mu(i)\)的板子
#include <bits/stdc++.h>
#define MAXN 3000000
#define LL long long
using namespace std;
const int mod = 1000000007;
const int inv2 = 500000004;
const int inv4 = 250000002;
const int inv6 = 166666668;
LL qp(LL a, LL b)
{
LL ret = 1;
while (b)
{
if (b & 1)
ret = ret * a % mod;
a = a * a % mod;
b >>= 1;
}
return ret;
}
LL inv(LL a)
{
return qp(a, mod - 2);
}
int prime[MAXN + 8];
bool notprime[MAXN + 8];
LL phi[MAXN + 8], phisum[MAXN + 8];
LL mu[MAXN + 8], musum[MAXN + 8];
LL iphisum[MAXN + 8];
LL imusum[MAXN + 8];
int ind = 0;
void getprime()
{
phi[1] = 1;
mu[1] = 1;
for (int i = 2; i <= MAXN; ++i)
{
if(!notprime[i])
{
prime[++ind] = i;
phi[i] = (i - 1) % mod;
mu[i] = -1;
}
for (int j = 1; j <= ind && i * prime[j] <= MAXN; ++j)
{
notprime[i * prime[j]] = 1;
if(i % prime[j] == 0)
{
phi[i * prime[j]] = phi[i] * prime[j] % mod;
mu[i * prime[j]] = 0;
break;
}
phi[i * prime[j]] = phi[i] * (prime[j] - 1) % mod;
mu[i * prime[j]] = -mu[i];
}
}
for (int i = 1; i <= MAXN; ++i)
{
phisum[i] = (phisum[i - 1] + phi[i]) % mod;
musum[i] = (musum[i - 1] + mu[i] + mod) % mod;
iphisum[i] = (iphisum[i - 1] + i * phi[i] % mod) % mod;
imusum[i] = (imusum[i - 1] + i * mu[i] % mod + mod) % mod;
}
}
unordered_map<LL,LL> _phi, _mu, _iphi, _imu;
LL Calcphi(LL n)
{
if (n <= MAXN)
return phisum[n];
auto it = _phi.find(n);
if (it != _phi.end())
return it->second;
LL last;
LL ans = (n % mod) * (((n % mod) + 1 + mod) % mod) % mod * inv2 % mod;
for (LL i = 2; i <= n; i = last + 1)
{
last = (n / (n / i));
ans = (ans - (((last % mod) - (i % mod) + 1 + mod) % mod * (Calcphi(n / i) % mod) % mod) + mod) % mod;
}
return _phi[n] = ans;
}
LL Calciphi(LL n)
{
if (n <= MAXN)
return iphisum[n];
auto it = _iphi.find(n);
if (it != _iphi.end())
return it->second;
LL ret = n * (n + 1) % mod * (2 * n + 1) % mod * inv6 % mod;
LL last;
for (LL i = 2; i <= n; i = last + 1)
{
last = n / (n / i);
ret = (ret - ((last + i) % mod * (last - i + 1) % mod * inv2 % mod + mod) % mod * Calciphi(n / i) % mod + mod) % mod;
}
return _iphi[n] = ret;
}
LL Calcmu(LL n)
{
if (n <= MAXN)
return musum[n];
auto it = _mu.find(n);
if (it != _mu.end())
return it->second;
LL ans = 1;
LL last;
for (LL i = 2; i <= n; i = last + 1)
{
last = n / (n / i);
ans = (ans - (last - i + 1) * Calcmu(n / i) % mod + mod) % mod;
}
return _mu[n] = ans;
}
LL Calcimu(int n)
{
if (n <= MAXN)
return imusum[n];
auto it = _imu.find(n);
if (it != _imu.end())
return it->second;
int last;
LL ret = 1;
for (int i = 2; i <= n; i = last + 1)
{
last = n / (n / i);
ret = (ret - ((1LL * last + i) % mod * (last - i + 1) % mod * inv2 % mod * Calcimu(n / i) % mod) + mod) % mod;
}
return _imu[n] = ret;
}
洲阁筛(还没学)
简单的写了个求1~n质数个数 -O2 N=1e11要2.5s
#pragma GCC optimize(2)
#include <bits/stdc++.h>
using namespace std;
#define LL long long
int L0[1000005],L1[1000005],sn,m,P[1000005],x,ga[1000005],gb[1000005];
LL n,a[1000005],b[1000005];
int main(){
n=100000000000;
// scanf("%lld",&n);
sn=sqrt(n);
for (int i=2;i<=sn;++i){
if (!P[i]) P[++m]=i;
for (int j=1;j<=m;++j){
x=P[j]*i;
if (x>sn) break;
P[x]=1;
if (i%P[j]==0) break;
}
}
for (int i=1;i<=sn;++i) a[i]=i,b[i]=n/i;
for (int i=1,j=1;i<=sn;++i){
while (j<=m&&P[j]*P[j]<=i) ++j;
L0[i]=j;
}
for (int i=sn,j=L0[sn];i;--i){
while (j<=m&&P[j]*P[j]<=b[i]) ++j;
L1[i]=j;
}
for (int i=1;i<=m;++i){
for (int j=1;j<=sn&&i<L1[j];++j){
LL y=j*P[i]; gb[j]=i;
if (y<=sn) b[j]-=b[y]+gb[y]+1-i;
else b[j]-=a[n/y]+ga[n/y]+1-i;
}
for (int j=sn;j&&i<L0[j];--j){
LL y=j/P[i]; ga[j]=i;
a[j]-=a[y]+ga[y]+1-i;
}
}
printf("%lld\n",b[1]+m-1);
return 0;
}
min25筛(学了一半)
待咕
类欧几里得算法
可快速计算:
\(f(a,b,c,n)=\sum\limits_{i=0}^n\lfloor\frac{ai+b}{c}\rfloor\)
\(g(a,b,c,n)=\sum\limits_{i=0}^ni\lfloor\frac{ai+b}{c}\rfloor\)
\(h(a,b,c,n)=\sum\limits_{i=0}^n\lfloor\frac{ai+b}{c}\rfloor^2\)
常用推导公式
\(a\leq \lfloor\frac{b}{c}\rfloor⇔ac\leq b\)
\(a< \lceil\frac{b}{c}\rceil⇔ac< b\)
\(a\geq \lceil\frac{b}{c}\rceil⇔ac\geq b\)
\(a> \lfloor\frac{b}{c}\rfloor⇔ac> b\)
\(\lfloor\frac{b}{c}\rfloor=\lfloor\frac{b+c-1}{c}\rfloor=\lceil\frac{b-c+1}{c}\rceil\)
以f的推导过程为例
1. \(a\geq c\)或\(b \geq c\)
\(\lfloor\frac{ai+b}{c}\rfloor=\lfloor\frac{(a\%c)i+b\%c}{c}\rfloor+\lfloor\frac{a}{c}\rfloor i+\lfloor\frac{b}{c}\rfloor\)
因此\(f(a,b,c,n)=f(a\%c,b\%c,c,n)+\frac{n(n+1)}{2}\lfloor\frac{a}{c}\rfloor+(n+1)\lfloor\frac{b}{c}\rfloor\)
2.\(a<c\)且\(b<c\)
\(a=0:f(a,b,c,n)=0\)
否则:
\(或g(a,b,c,n)=g(a\%c,b\%c,c,n)+\frac{n(n+1)(2n+1)}{6}\lfloor\frac{a}{c}\rfloor+\frac{n(n+1)}{2}\lfloor\frac{b}{c}\rfloor\ (a\geq c或b\geq c)\)
\(否则g(a,b,c,n)=\frac{\lfloor\frac{an+b}{c}\rfloor n(n+1)-f(c,-b+c-1,a,\lfloor\frac{an+b}{c}\rfloor-1)-h(c,-b+c-1,a,\lfloor\frac{an+b}{c}\rfloor-1)}{2}(否则)\)
\(或h(a,b,c,n)=h(a\%c,b\%c,c,n)+\frac{n(n+1)(2n+1)}{6}\lfloor\frac{a}{c}\rfloor^2+(n+1)\lfloor\frac{b}{c}\rfloor^2+2\lfloor\frac{b}{c}\rfloor f(a\%c,b\%c,c,n)\\+2\lfloor\frac{a}{c}\rfloor g(a\%c,b\%c,c,n)+\lfloor\frac{a}{c}\rfloor\lfloor\frac{b}{c}\rfloor n(n+1)(a\geq c或b\geq c)\)
\(否则h(a,b,c,n)=\lfloor\frac{an+b}{c}\rfloor(\lfloor\frac{an+b}{c}\rfloor+1)n-2g(c,-b+c-1,a,\lfloor\frac{an+b}{c}\rfloor-1)\\-2f(c,-b+c-1,a,\lfloor\frac{an+b}{c}\rfloor-1)-f(a,b,c,n)(否则)\)
/*洛谷P5170,输入T组,一开始写了个递归暴力,给我T飞掉了,我人都傻了*/
#include <bits/stdc++.h>
#define LL long long
#define MAXN 2000
#define x first
#define y second
#define pii pair<long long, long long>
using namespace std;
const LL mod = 998244353;
const LL inv2 = 499122177;
const LL inv6 = 166374059;
const LL inv4 = 748683265;
struct Node
{
LL f, g, h;
Node(){}
Node(LL _f, LL _g, LL _h)
{
f = _f; g = _g; h = _h;
}
};
Node solve(LL a, LL b, LL c, LL n)
{
if (a == 0)
{
n %= mod;
LL bc = (b/c) % mod;
return Node((n+1)*bc%mod, n*(n+1)%mod*inv2%mod*bc%mod, (n+1)*bc%mod*bc%mod);
}
LL t = (a*n+b)/c;
LL ac=(a/c)%mod, bc=(b/c)%mod;
if (a >= c || b >= c)
{
Node n1 = solve(a%c,b%c,c,n);
return Node(( n1.f + n*(n+1)%mod*inv2%mod*ac%mod + (n+1)*bc%mod) % mod,
( n1.g + n*(n+1)%mod*(2*n+1)%mod*inv6%mod*ac%mod % mod + n*(n+1)%mod*inv2%mod*bc%mod) % mod,
( n1.h + n*(n+1)%mod*(2*n+1)%mod*inv6%mod*ac%mod*ac%mod + (n+1)*bc%mod*bc%mod + 2*bc%mod*n1.f%mod + 2*ac%mod*n1.g%mod + ac*bc%mod*n%mod*(n+1)%mod) % mod
);
}
else
{
Node n2 = solve(c,c-b-1,a,t-1);
t %= mod;
return Node((n * t % mod - n2.f + mod) % mod,
inv2 * ((t*n%mod*(n+1)%mod - n2.f - n2.h + 2*mod) % mod) % mod,
(t*(t+1)%mod*n%mod - 2*n2.g%mod - 2*n2.f%mod - (n * t % mod - n2.f + mod) % mod + 4*mod) % mod
);
}
}
int main()
{
LL a, b, c, n;
int T;
scanf("%d", &T);
while (T--)
{
scanf("%lld %lld %lld %lld", &n, &a, &b, &c);
Node ans = solve(a, b, c, n);
printf("%lld %lld %lld\n", ans.f, ans.h, ans.g);
}
return 0;
}
拉格朗日插值法
朴素算法
给定\(k+1\)个点\((x_0,y_0),...,(x_k,y_k)\)表示的\(k\)次多项式\(L(x)\),给出\(a\)求\(L(a)\).
记\(\ell_j(x)=\prod\limits_{i=0,i\neq j}^k\frac{x-x_i}{x_j-x_i}\),则有\(L(x)=\sum\limits_{j=0}^ky_j\ell_j(x)\). \(O(n^2)\).
/*洛谷P4781 给出n,X,接下来给出n对(x,y),求L(X).*/
#include <bits/stdc++.h>
#define LL long long
#define MAXN 2000
#define x first
#define y second
#define pii pair<long long, long long>
using namespace std;
const int mod = 998244353;
int n;
LL X;
pii arr[MAXN + 8];
LL qp(LL a, LL b)
{
LL ret = 1;
while (b)
{
if (b & 1)
ret = ret * a % mod;
a = a * a % mod;
b >>= 1;
}
return ret;
}
LL inv(LL a)
{
return qp(a, mod - 2);
}
LL calcell(int j, LL X)
{
LL ret = 1;
for (int i = 1; i <= n; ++i)
{
if (i == j || arr[i].x == arr[j].x)
continue;
ret = ret * (X - arr[i].x + mod) % mod * inv((arr[j].x - arr[i].x + mod) % mod) % mod;
}
return ret;
}
LL calcF(LL X)
{
LL ret = 0;
for (int i = 1; i <= n; ++i)
{
ret = (ret + arr[i].y * calcell(i, X) % mod) % mod;
}
return ret;
}
int main()
{
cin >> n >> X;
for (int i = 1; i <= n; ++i)
{
cin >> arr[i].x >> arr[i].y;
}
cout << calcF(X) << endl;
return 0;
}
重心拉格朗日插值法
插入新点时单次插入优化至\(O(n)\):
记\(\ell(x)=\prod\limits_{i=0}^k(x-x_i)\),\(w_j=\frac{1}{\prod_{i=0,i\neq j}^k(x_j-x_i)}\),则\(\ell_j(x)=\ell(x)\frac{w_j}{x-x_j}\).
则\(L(x)=\ell(x)\sum\limits_{j=0}^k\frac{w_j}{x-x_j}y_j\).
对多项式\(g(x)=1\)插值,有\(\forall x,g(x)=\ell(x)\sum\limits_{j=0}^k\frac{w_j}{x-x_j}\),因此\(L(x)\)两边同除以\(g(x)\),有
\(L(x)=\frac{\sum_{j=0}^k\frac{w_j}{x-x_j}y_j}{\sum_{j=0}^k\frac{w_j}{x-x_j}}\).插入新值\(x_{k+1}\)时,只需对\(w_j\)乘上\(x_j-x_{k+1}\)即可更新\(w\),更新所有\(w\)总共需要\(O(n)\).
/*依然是上面那个模板题*/
#include <bits/stdc++.h>
#define LL long long
#define MAXN 2000
#define x first
#define y second
#define pii pair<long long, long long>
using namespace std;
const int mod = 998244353;
int n;
LL X;
pii arr[MAXN + 8];
LL w[MAXN + 8];
LL qp(LL a, LL b)
{
LL ret = 1;
while (b)
{
if (b & 1)
ret = ret * a % mod;
a = a * a % mod;
b >>= 1;
}
return ret;
}
LL inv(LL a)
{
return qp(a, mod - 2);
}
void insert(int i)
{
w[i] = 1;
for (int j = 1; j < i; ++j)
{
if (arr[i].x == arr[j].x)
continue;
w[i] = w[i] * (arr[i].x - arr[j].x + mod) % mod;
w[j] = w[j] * inv(arr[j].x - arr[i].x + mod) % mod;
}
w[i] = inv(w[i]);
return;
}
LL query(int n)
{
LL fz = 0, fm = 0;
for (int i = 1; i <= n; ++i)
{
fz = (fz + w[i] * inv((X - arr[i].x + mod) % mod) % mod * arr[i].y % mod) % mod;
fm = (fm + w[i] * inv((X - arr[i].x + mod) % mod) % mod) % mod;
}
return fz * inv(fm) % mod;
}
int main()
{
cin >> n >> X;
for (int i = 1; i <= n; ++i)
{
cin >> arr[i].x >> arr[i].y;
insert(i);
}
cout << query(n) << endl;
return 0;
}
佩尔方程
\(x^2-Dy^2=1,D\)不是平方数。
设\((x_1,y_1)\)是使\(x_1\)最小的解,则通解为
\(x_k+\sqrt{D}y_k=(x_1+\sqrt{D}y_1)^k\).
\(x^2-Dy^2=1,D\)不是平方数。
设\((x_0,y_0)\)是使\(x_0\)最小的解,则通解为
\(x_k+\sqrt{D}y_k=(x_0+\sqrt{D}y_0)^{2k+1}\).
佩尔方程(\(x^2-Dy^2=1\))的最小正整数解
bool PQA(LL D, LL &p, LL &q) {//来自于PQA算法的一个特例
const double eps = 1e-6;
LL d = sqrt(D);
if (abs(d * d - D) <= eps)
return false;//这里是判断佩尔方程有没有解
LL u = 0, v = 1, a = int(sqrt(D)), a0 = a, lastp = 1, lastq = 0;
p = a, q = 1;
do {
u = a * v - u;
v = (D - u * u) / v;
a = (a0 + u) / v;
LL thisp = p, thisq = q;
p = a * p + lastp;
q = a * q + lastq;
lastp = thisp;
lastq = thisq;
} while ((v != 1 && a <= a0));//这里一定要用do~while循环
p = lastp;
q = lastq;
//这样求出后的(p,q)是p2 – D * q2 = (-1)k的解,也就是说p2 – D * q2可能等于1也可能等于-1,如果等于1,(p,q)就是解,如果等于-1还要通过(p2 + D * q2,2 * p * q)来求解,如下
if (p * p - D * q * q == -1) {
p = lastp * lastp + D * lastq * lastq;
q = 2 * lastp * lastq;
}
return true;
}
正整数的拆分方案数\(O(n\sqrt{n})\)
void init() {
for(int i = 0; i < 1000; ++i) {
g[i << 1] = 1LL * i * (i * 3 - 1) / 2;
g[i << 1 | 1] = 1LL * i * (i * 3 + 1) / 2;
}
f[0] = 1, f[1] = 1, f[2] = 2;
for(int i = 3; i < N; ++i) {
f[i] = 0;
int k = 0;
for(int j = 2; g[j] <= i; ++j) {
if(k & 2)
f[i] = ((f[i] - f[i - g[j]]) % mod + mod) % mod;
else
f[i] = ((f[i] + f[i - g[j]]) % mod + mod) % mod;
k++; k %= 8;
}
}
}
/*g是五边形数,f是拆分方案数*/
/*若限制每个数大于等于3,则 ans=f[i]+f[i-3]-f[i-1]-f[i-2]*/
博弈论
常用简单博弈模型
巴什博弈
一堆n个石子,两人轮流取,每次最多取m个,取完最后一颗的人胜。
解:\(n\%(m+1) = 0\)时,先手必败,否则必胜。
斐波那契博弈
一堆n个石子,先手取任意多个,但不能取完,之后每次取的石子数只能在1到对手刚取的石子数的2倍之间(闭区间),取完最后一颗的人胜。
解:\(n\)是斐波那契数时,先手必败,否则必胜。
威佐夫博弈
两堆各若干个石子,每次可以一堆取任意多个,或两堆取同样多个,取完最后一颗的胜。
解:\((0,0),(1,2),(3,5),(4,7),(6,10)....\)先手必败。一般地,对奇异局势(先手必败局势)\((a,b)(a\le b),\)有$a=\lfloor(b-a)\frac{\sqrt 5 + 1}{2}\rfloor $
或者说,第\(k\)组奇异局势\((a_k,b_k)\)有\(a_k=\lfloor\frac{\sqrt 5 + 1}{2}k\rfloor, b_k=a_k +k\).
尼姆博弈
若干堆石子,每堆若干个,每次选一堆取至少一个,取完最后一颗的胜。
解:\(\bigoplus_{i=1}^n a[i] = 0\)时,先手必败,否则必胜。
阶梯博弈
若干堆石子,堆标号依次\(1,2,3,...,n\)每次选一堆(标号为\(i\))的若干个放到标号\(i-1\)的堆(可以放到标号为0的堆,称作地面),把最后一颗石子放到\(0\)号堆的胜。
解:等价于保留奇数号,忽视偶数号的Nim博弈。