求积性前缀和,即
S ( n ) = ∑ i = 1 n f ( i ) S(n)=\sum_{i=1}^nf(i) S ( n ) = ∑ i = 1 n f ( i ) .
某些不是积性的函数也可以,只要能找到一个合适的
g g g 。
对于任意两个数论函数
f , g f,g f , g ,有
∑ i = 1 n ( f ∗ g ) ( i ) = ∑ i = 1 n g ( i ) S ( ⌊ n i ⌋ ) \sum_{i=1}^n(f*g)(i)=\sum_{i=1}^ng(i)S(\lfloor\frac{n}{i}\rfloor) ∑ i = 1 n ( f ∗ g ) ( i ) = ∑ i = 1 n g ( i ) S (⌊ i n ⌋) .
∑ i = 1 n ( f ∗ g ) ( i ) = ∑ i = 1 n ∑ d ∣ n g ( d ) f ( i d ) = ∑ i = 1 n g ( i ) ∑ j = 1 ⌊ n / i ⌋ f ( j ) = ∑ i = 1 n g ( i ) S ( ⌊ n i ⌋ ) \sum_{i=1}^n(f*g)(i)=\sum_{i=1}^n\sum_{d|n}g(d)f(\frac{i}{d})=\sum_{i=1}^ng(i)\sum_{j=1}^{\lfloor n/i\rfloor}f(j)=\sum_{i=1}^ng(i)S(\lfloor\frac{n}{i}\rfloor) ∑ i = 1 n ( f ∗ g ) ( i ) = ∑ i = 1 n ∑ d ∣ n g ( d ) f ( d i ) = ∑ i = 1 n g ( i ) ∑ j = 1 ⌊ n / i ⌋ f ( j ) = ∑ i = 1 n g ( i ) S (⌊ i n ⌋)
如果可以找到一个
g ( n ) g(n) g ( n ) 可以快速地求出
f ∗ g f*g f ∗ g 的前缀和及
g g g 的前缀和,则
g ( 1 ) S ( n ) = ∑ i = 1 n ( f ∗ g ) ( i ) − ∑ i = 2 n g ( i ) S ( ⌊ n i ⌋ ) g(1)S(n)=\sum_{i=1}^n(f*g)(i)-\sum_{i=2}^ng(i)S(\lfloor\frac{n}{i}\rfloor) g ( 1 ) S ( n ) = ∑ i = 1 n ( f ∗ g ) ( i ) − ∑ i = 2 n g ( i ) S (⌊ i n ⌋) 。
后半部分可以递归地上数论分块,求得的若干
S ( n ) S(n) S ( n ) 需要记忆化来保证复杂度。
总复杂度
O ( n 3 / 4 ) O(n^{3/4}) O ( n 3/4 ) ,若预处理前
O ( n 2 / 3 ) O(n^{2/3}) O ( n 2/3 ) 项,则复杂度降至
O ( n 2 / 3 ) O(n^{2/3}) O ( n 2/3 ) 。
多项式乘法
大家伙。
计算两个次数为
n n n 的多项式
F ( x ) , G ( x ) F(x),G(x) F ( x ) , G ( x ) 的乘积。
系数表示法计算乘积是
O ( n 2 ) O(n^2) O ( n 2 ) 的,而点值表示法仅需要
O ( n ) O(n) O ( n ) 。尝试把两个多项式转化为点值表示,相乘再转化回去。
由拉格朗日插值法,
n n n 次多项式需要
n + 1 n+1 n + 1 个点来唯一确定。在低于
O ( n 2 ) O(n^2) O ( n 2 ) 的复杂度下对于
n + 1 n+1 n + 1 个点计算任意一个多项式
F ( x ) F(x) F ( x ) 看起来也很扯,只能寄希望于这些点满足一些奇妙性质。
以下认为
n + 1 = 2 k n+1=2^k n + 1 = 2 k ,如果不是就用 0 补齐高次项。
考虑将
n n n 次单位根
1 , ω n , ω n 2 , ⋯ , ω n n − 1 1,\omega_n,\omega_n^2,\cdots,\omega_n^{n-1} 1 , ω n , ω n 2 , ⋯ , ω n n − 1 代入
F ( x ) F(x) F ( x ) 求值。
将
F ( x ) F(x) F ( x ) 的系数奇偶分类,奇数分给
A A A ,偶数分给
B B B 。则
F ( x ) = A ( x 2 ) + x B ( x 2 ) F(x)=A(x^2)+xB(x^2) F ( x ) = A ( x 2 ) + x B ( x 2 ) 。注意到:
F ( ω n k ) = A ( ( ω n k ) 2 ) + ω n k B ( ( ω n k ) 2 ) = A ( ω n / 2 k ) + ω n k B ( ω n / 2 k ) F(\omega_n^k)=A((\omega_n^k)^2)+\omega_n^kB((\omega_n^k)^2)=A(\omega_{n/2}^k)+\omega_n^kB(\omega_{n/2}^k) F ( ω n k ) = A (( ω n k ) 2 ) + ω n k B (( ω n k ) 2 ) = A ( ω n /2 k ) + ω n k B ( ω n /2 k )
因为
ω n k + n / 2 = − ω n k \omega_{n}^{k+n/2}=-\omega_{n}^{k} ω n k + n /2 = − ω n k ,所以
F ( ω n k + n / 2 ) = A ( ω n / 2 k ) − ω n k B ( ω n / 2 k ) F(\omega_{n}^{k+n/2})=A(\omega_{n/2}^k)-\omega_n^kB(\omega_{n/2}^k) F ( ω n k + n /2 ) = A ( ω n /2 k ) − ω n k B ( ω n /2 k )
所以如果分别获得了
A ( ω n / 2 k ) , B ( ω n / 2 k ) A(\omega_{n/2}^k),B(\omega_{n/2}^k) A ( ω n /2 k ) , B ( ω n /2 k ) 的值,我们就可以一起算出
F ( ω n k ) F(\omega_n^k) F ( ω n k ) 和
F ( ω n k + n / 2 ) F(\omega_{n}^{k+n/2}) F ( ω n k + n /2 ) 。而计算
A ( ω n / 2 k ) , B ( ω n / 2 k ) A(\omega_{n/2}^k),B(\omega_{n/2}^k) A ( ω n /2 k ) , B ( ω n /2 k ) 是一个规模更小的子问题,递归求解是可以的。
这是
O ( n log n ) O(n\log n) O ( n log n ) 的。将系数表示法转化为点值表示法仅需要将所有单位根取共轭(虚部取反),进行相同的操作并将结果的每一项
× 1 n \times \frac{1}{n} × n 1 。
鬼知道为什么。
浅证一下吧。设
{ a i } \{a_i\} { a i } 为
F F F 各项系数,
b i = F ( ω n i ) = ∑ k = 0 n − 1 a k ( ω n i ) k b_i=F(\omega_n^i)=\sum_{k=0}^{n-1}a_k(\omega_n^i)^k b i = F ( ω n i ) = ∑ k = 0 n − 1 a k ( ω n i ) k ,
c i = 1 n ∑ k = 0 n − 1 b k ( ω n − i ) k c_i=\frac{1}{n}\sum_{k=0}^{n-1}b_k(\omega_n^{-i})^k c i = n 1 ∑ k = 0 n − 1 b k ( ω n − i ) k 。那我们要证
c i = a i c_i=a_i c i = a i 。直接展开。
c i = 1 n ∑ k = 0 n − 1 ∑ j = 0 n − 1 a j ( ω n k ) j ( ω n − i ) k = 1 n ∑ j = 0 n − 1 a j ∑ k = 0 n − 1 ( ω n j − i ) k = 1 n ∑ j = 0 n − 1 [ j = i ] n ⋅ a j = a i \begin{aligned}
c_i&=\frac{1}{n}\sum_{k=0}^{n-1}\sum_{j=0}^{n-1}a_j(\omega_n^k)^j(\omega_n^{-i})^k
\\
&=\frac{1}{n}\sum_{j=0}^{n-1}a_j\boxed{\sum_{k=0}^{n-1}(\omega_n^{j-i})^k}
\\
&=\frac{1}{n}\sum_{j=0}^{n-1}\boxed{[j=i]n}\cdot a_j
\\
&=a_i
\end{aligned} c i = n 1 k = 0 ∑ n − 1 j = 0 ∑ n − 1 a j ( ω n k ) j ( ω n − i ) k = n 1 j = 0 ∑ n − 1 a j k = 0 ∑ n − 1 ( ω n j − i ) k = n 1 j = 0 ∑ n − 1 [ j = i ] n ⋅ a j = a i
关于单位根求和,既可以直接套等比求和,也可以在单位圆上看一眼。都行,无所谓。
具体实现上,会先找到
F ( i ) F(i) F ( i ) 经过一串区间内奇偶分类的操作后会落在哪里。将初始位置二进制反转后就可以得到最终位置。鬼知道为什么。
源码大概长这个样子。
CPP int rev[N];
void shuffle (cmplx a[],int len) {
for (int i=1 ;i<len;i++){
rev[i]=rev[i>>1 ]>>1 |(i&1 ?(len>>1 ):0 );
}
for (int i=0 ;i<len;i++) if (i<rev[i]) swap (a[i],a[rev[i]]);
}
void fast_fast_TLE (cmplx a[],int len,int flag) {
shuffle (a,len);
for (int k=2 ;k<=len;k<<=1 ){
int wid=k>>1 ;
cmplx w1={cos (2 *pi/k),flag*sin (2 *pi/k)};
for (int i=0 ;i<len;i+=k){
cmplx w={1 ,0 };
for (int j=i;j<i+wid;j++){
cmplx p=a[j],q=a[j+wid]*w;
a[j]=p+q; a[j+wid]=p-q;
w*=w1;
}
}
}
if (flag==-1 ) for (int i=0 ;i<len;i++) a[i].real/=len;
}
NTT
考虑我们用到了上述单位根的什么性质。
ω 2 n 2 k = ω n k \omega_{2n}^{2k}=\omega_{n}^k ω 2 n 2 k = ω n k
ω 2 n k + n = − ω 2 n k \omega_{2n}^{k+n}=-\omega_{2n}^k ω 2 n k + n = − ω 2 n k
( ω n k ) 2 = ω n 2 k (\omega_{n}^k)^2=\omega_{n}^{2k} ( ω n k ) 2 = ω n 2 k
我们发现如果设
ω n k = g ( P − 1 ) k n \omega_{n}^k=g^{(P-1)\frac{k}{n}} ω n k = g ( P − 1 ) n k ,在
m o d P \bmod P mod P 下好像也满足上述性质。
(一个数
p p p 的原根
g g g 满足
g k g^{k} g k 是纯循环(由欧拉定理,
gcd ( g , p ) = 1 \gcd(g,p)=1 g cd( g , p ) = 1 )且循环节为
φ ( k ) \varphi(k) φ ( k ) 。
p p p 是质数时即代表
g 0 , g 1 , ⋯ , g P − 2 g^0,g^1,\cdots,g^{P-2} g 0 , g 1 , ⋯ , g P − 2 遍历
[ 1 , P ) [1,P) [ 1 , P ) 。)
所以我们将上面的
ω \omega ω 全文替换便得到了单位根。
注意并不是所有模数都可以
NTT \operatorname{NTT} NTT 。求单位根的时候我们将
( P − 1 ) (P-1) ( P − 1 ) 除以了
n n n ,我们需要让其是整数。模数是
998244353 = 119 × 2 23 + 1 998244353=119\times 2^{23}+1 998244353 = 119 × 2 23 + 1 ,允许我们递归
23 23 23 层
NTT \operatorname{NTT} NTT ,即多项式最高不能超过
2 23 2^{23} 2 23 次。实际使用时不用担心,因为次数高了最先爆炸的是时限而不是正确性。
写出来大概长这个样子:
CPP void faster_faster_TLE (int a[],int len,int flag) {
shuffle (a,len);
for (int k=2 ;k<=len;k<<=1 ){
int wid=k>>1 ;
int w1=Pow (G,(flag*(P-1 )/k+(P-1 ))%(P-1 ));
for (int i=0 ;i<len;i+=k){
int w=1 ;
for (int j=i;j<i+wid;j++){
int p=a[j],q=a[j+wid]*w%P;
a[j]=(p+q)%P; a[j+wid]=(p-q+P)%P;
w*=w1; w%=P;
}
}
}
if (flag==-1 ){
int Ilen=Inv (len);
for (int i=0 ;i<len;i++) a[i]*=Ilen,a[i]%=P;
}
}
停时
记一个随机现象中所有可能发生的结果(称为样本点)组成的集合为
Ω \Omega Ω ,记
ω \omega ω 是一个随机试验的结果:
ω \omega ω 是确定取一个
∈ Ω \in\Omega ∈ Ω 中的元素,但是在正式试验之前我们不知道它具体取什么。
记
X ( ω ) X(\omega) X ( ω ) 为一个随机变量。
X X X 本质上是一个函数,将
样本点 映射到
数值 上。这样来说
ω \omega ω 就很像一个生成种子。
记一个随机过程为一系列的随机变量的集合(序列?),每个元素对应着一个时刻
t t t 。比如说像这样:
{ X 0 ( ω ) , X 1 ( ω ) , X 2 ( ω ) , ⋯ } \{X_0(\omega),X_1(\omega),X_2(\omega),\cdots\} { X 0 ( ω ) , X 1 ( ω ) , X 2 ( ω ) , ⋯ } 。注意所有
X X X 的取值均依赖于一个
ω \omega ω ,即一旦
ω \omega ω 确定了整个随机过程也定下来了。把
ω \omega ω 略掉不写是可以的。
定义鞅为“未来的值的期望等于当前的值”的随机过程。形式化来讲,若随机过程
{ X t } t ≥ 0 \{X_t\}_{t\geq 0} { X t } t ≥ 0 对于每一个时刻
t t t 满足
E ( X t + 1 ∣ F t ) = X t \mathbb{E}(X_{t+1}|\mathcal{F}_t)=X_t E ( X t + 1 ∣ F t ) = X t 。其中
F t \mathcal F_t F t 表示截至
t t t 时刻所获得的随机过程的信息(?),即
X 0 , X 1 , ⋯ , X t X_0,X_1,\cdots,X_t X 0 , X 1 , ⋯ , X t 。(应该是,只确定了这些
X X X 的值,而
ω \omega ω 仍然不确定?这句话是自己瞎猜的。)
一个鞅的停时
τ \tau τ 是一个随机变量,满足
∀ t ≥ 0 , { τ = t } ∈ F t \forall t\geq 0,\{\tau=t\}\in\mathcal F_t ∀ t ≥ 0 , { τ = t } ∈ F t ,即到任意时刻时我都知道当前该不该停下来。
停时定理
叙述如下:
设
{ X t } t ≥ 0 \{X_t\}_{t\geq 0} { X t } t ≥ 0 是一个鞅,
τ \tau τ 是停时,如果满足以下条件
之一 :
【censored】
【censored】
【censored】
则
E [ X τ ] = E [ X 0 ] \mathbb{E}[X_\tau]=\mathbb E[X_0] E [ X τ ] = E [ X 0 ] 。
(条件是什么不重要,即使知道也不一定会判断是否满足,直接当成永远成立吧。)
势能函数
目的是求解一个随机过程(应该不需要是鞅)的停时的期望
E ( τ ) \mathbb E(\tau) E ( τ ) 。
构造一个势能函数
Φ ( X ) \Phi(X) Φ ( X ) ,满足:
E ( Φ ( X t + 1 ) − Φ ( X t ) ∣ F n ) = − 1 \mathbb E(\Phi(X_{t+1})-\Phi(X_t)|\mathcal F_n)=-1 E ( Φ ( X t + 1 ) − Φ ( X t ) ∣ F n ) = − 1 ,
∀ X τ \forall X_\tau ∀ X τ ,Φ ( X τ ) \Phi(X_\tau) Φ ( X τ ) 均相等,
Φ ( X i ) = Φ ( X τ ) \Phi(X_i)=\Phi(X_\tau) Φ ( X i ) = Φ ( X τ ) 当且仅当 i = τ i=\tau i = τ 。
人话就是势能函数每次期望
− 1 -1 − 1 ,停止时到达一个常数。
(其实没有理解为什么有第三条?防止
i i i 被误认为是停时?但是不知道为什么这样被误认会有问题(
然后由于期望的线性性,裂项就(应该?)得到了
E ( τ ) = E ( Φ ( X 0 ) ) − Φ ( X τ ) \mathbb E(\tau)=\mathbb E(\Phi(X_0))-\Phi(X_\tau) E ( τ ) = E ( Φ ( X 0 )) − Φ ( X τ ) 。停时的期望就转化成了势能函数的期望。所以我们需要做的就是确定一个美好的势能函数然后竭尽全力计算其初态的期望。
以下设
s s s 为球的总数,
c i c_i c i 为颜色为
i i i 的球的数量。
对于当前这道题,
S S S 表示球的颜色状态,如果让势能函数
Φ ( S ) = ∑ f ( c i ) \Phi(S)=\sum f(c_i) Φ ( S ) = ∑ f ( c i ) ,其中
f f f 是某个函数,
c i c_i c i 是颜色为
i i i 的球的个数,则终止状态的势能函数值自然全相同。
因为期望的可加性,可直接考虑单个
f ( c ) f(c) f ( c ) 转移一次得到的东西的期望:
f ( c ) → c ( s − c ) s ( s − 1 ) ( f ( c − 1 ) + f ( c + 1 ) ) + ( 1 − 2 c ( s − c ) s ( s − 1 ) f ( c ) ) f(c)\rightarrow\frac{c(s-c)}{s(s-1)}(f(c-1)+f(c+1))+\left(1-\frac{2c(s-c)}{s(s-1)}f(c)\right) f ( c ) → s ( s − 1 ) c ( s − c ) ( f ( c − 1 ) + f ( c + 1 )) + ( 1 − s ( s − 1 ) 2 c ( s − c ) f ( c ) )
分别讨论了用当前颜色覆盖另一颜色的球,用另一颜色覆盖当前颜色的球,两个球都是/都不是当前颜色的情况。
设箭头右边的一坨为
f ′ ( c ) f'(c) f ′ ( c ) 。然后我们需要让
F ( S ) F(S) F ( S ) 转移一次期望
− 1 -1 − 1 ,所以即对于任意
∑ c i = s \sum c_i=s ∑ c i = s 的
{ c } \{c\} { c } ,都有
∑ ( f ( c i ) − f ′ ( c i ) ) = 1 \sum(f(c_i)-f'(c_i))=1 ∑ ( f ( c i ) − f ′ ( c i )) = 1 。
f ′ ( c ) f'(c) f ′ ( c ) 只能为
f ( c ) − c s f(c)-\frac{c}{s} f ( c ) − s c 。
回到上面超长的柿子,把等式套进去解,得到
2 f ( c ) = f ( c + 1 ) + f ( c − 1 ) + s − 1 s − c 2f(c)=f(c+1)+f(c-1)+\frac{s-1}{s-c} 2 f ( c ) = f ( c + 1 ) + f ( c − 1 ) + s − c s − 1
设
d ( c ) = f ( c + 1 ) − f ( c ) d(c)=f(c+1)-f(c) d ( c ) = f ( c + 1 ) − f ( c ) ,则改写成
d ( c ) = d ( c − 1 ) − s − 1 s − c d(c)=d(c-1)-\frac{s-1}{s-c} d ( c ) = d ( c − 1 ) − s − c s − 1 。
规定
f ( 0 ) = 0 , d ( 0 ) = − 1 f(0)=0,d(0)=-1 f ( 0 ) = 0 , d ( 0 ) = − 1 ,然后快乐地递推就可以得到每一个
f ( c i ) f(c_i) f ( c i ) 了。现在还需要求
f ( s ) f(s) f ( s ) 。
f ( s ) = ∑ i = 0 s − 1 d ( i ) = ∑ i = 0 s − 1 d ( 0 ) + ∑ j = 1 i s − 1 s − j = ( s − 1 ) 2 + s × d ( 0 ) f(s)=\sum_{i=0}^{s-1}d(i)=\sum_{i=0}^{s-1}d(0)+\sum_{j=1}^i\frac{s-1}{s-j}=(s-1)^2+s\times d(0) f ( s ) = ∑ i = 0 s − 1 d ( i ) = ∑ i = 0 s − 1 d ( 0 ) + ∑ j = 1 i s − j s − 1 = ( s − 1 ) 2 + s × d ( 0 )
做完辣。
还是记
Φ ( S ) = ∑ f ( s i ) \Phi(S)=\sum f(s_i) Φ ( S ) = ∑ f ( s i ) 。设一次操作选择了
x x x 个颜色
s i s_i s i 的球,
y y y 个颜色不是
s i s_i s i 的球。发生的概率为
( s i x ) ( n − s i y ) 2 n \frac{\binom{s_i}{x}\binom{n-s_i}{y}}{2^n} 2 n ( x s i ) ( y n − s i ) 。将一个被选中的球染成
s i s_i s i 的概率为
x + y n \frac{x+y}{n} n x + y 。于是我们可以得到下式:
− 1 + ∑ i = 1 n f ( s i ) = ∑ i = 1 n ∑ x = 0 s i ∑ y = 0 n − s i ( s i x ) ( n − s i y ) 2 n ⋅ ( x + y n f ( s i − x + 1 ) + n − x − y n f ( s i − x ) ) -1+\sum_{i=1}^nf(s_i)=\sum_{i=1}^n\sum_{x=0}^{s_i}\sum_{y=0}^{n-s_i}\frac{\binom{s_i}{x}\binom{n-s_i}{y}}{2^n}\cdot\left(\frac{x+y}{n}f(s_i-x+1)+\frac{n-x-y}{n}f(s_i-x)\right) − 1 + ∑ i = 1 n f ( s i ) = ∑ i = 1 n ∑ x = 0 s i ∑ y = 0 n − s i 2 n ( x s i ) ( y n − s i ) ⋅ ( n x + y f ( s i − x + 1 ) + n n − x − y f ( s i − x ) )
(我知道你很急但你先别急(
然后因为对于所有
s i s_i s i 的组合都应成立所以
i i i 的求和可以削去。(并把
− 1 -1 − 1 变成
∑ i = 1 n s i n \sum_{i=1}^n \frac{s_i}{n} ∑ i = 1 n n s i 的形式。)
省略下标
i i i ,直接记
s i = s s_i=s s i = s ,变换然后把混在一起的变量分开一下:
f ( s ) − s n = 1 n 2 n ∑ x = 0 s ( s x ) f ( s − x + 1 ) ∑ y = 0 n − s ( x + y ) ( n − s y ) + 1 n 2 n ∑ x = 0 s ( s x ) f ( s − x ) ∑ y = 0 n − s ( n − x − y ) ( n − s y ) \begin{aligned}
f(s)-\frac{s}{n}&=
\frac{1}{n2^n}\sum_{x=0}^s\binom{s}{x}f(s-x+1)\boxed{\sum_{y=0}^{n-s}(x+y)\binom{n-s}{y}}\\
&+\frac{1}{n2^n}\sum_{x=0}^s\binom{s}{x}f(s-x)\boxed{\sum_{y=0}^{n-s}(n-x-y)\binom{n-s}{y}}
\end{aligned} f ( s ) − n s = n 2 n 1 x = 0 ∑ s ( x s ) f ( s − x + 1 ) y = 0 ∑ n − s ( x + y ) ( y n − s ) + n 2 n 1 x = 0 ∑ s ( x s ) f ( s − x ) y = 0 ∑ n − s ( n − x − y ) ( y n − s )
然后你知道
∑ i = 0 n ( n i ) = 2 n , ∑ i = 0 n i ( n i ) = ( n − 1 ) 2 n − 1 \sum_{i=0}^n\binom ni=2^n,\sum_{i=0}^ni\binom ni=(n-1)2^{n-1} ∑ i = 0 n ( i n ) = 2 n , ∑ i = 0 n i ( i n ) = ( n − 1 ) 2 n − 1 ,于是框起来的部分就无了。
具体地,第一行的尾巴是
( n + 2 x − s ) 2 n − s − 1 (n+2x-s)2^{n-s-1} ( n + 2 x − s ) 2 n − s − 1 ,第二行的尾巴是
( n − 2 x + s ) 2 n − s − 1 (n-2x+s)2^{n-s-1} ( n − 2 x + s ) 2 n − s − 1 。
回代整理可得:
RHS = 1 n 2 s + 1 ∑ x = 0 s ( s x ) ( ( n + 2 x − s ) f ( s − x + 1 ) + ( n − 2 x + s ) f ( s − x ) ) \text{RHS}=\frac{1}{n2^{s+1}}\sum_{x=0}^s\binom{s}{x}\big((n+2x-s)f(s-x+1)+(n-2x+s)f(s-x)\big) RHS = n 2 s + 1 1 ∑ x = 0 s ( x s ) ( ( n + 2 x − s ) f ( s − x + 1 ) + ( n − 2 x + s ) f ( s − x ) )
然后你要把和式里面的
f ( s + 1 ) f(s+1) f ( s + 1 ) 单拎出来因为这才是你想求的:(顺便,
x ′ ← s − x x'\leftarrow s-x x ′ ← s − x )
− f ( s + 1 ) = n + s − n 2 s + 1 n − s f ( s ) + s n − s 2 s + 1 + 1 n − s ∑ x = 0 s − 1 ( s x ) ( ( n − 2 x + s ) f ( x + 1 ) + ( n + 2 x − s ) f ( x ) ) \begin{aligned}
-f(s+1)&=\frac{n+s-n2^{s+1}}{n-s}f(s)+\frac{s}{n-s}2^{s+1}\\
&+\frac{1}{n-s}\sum_{x=0}^{s-1}\binom sx\big((n-2x+s)f(x+1)+(n+2x-s)f(x)\big)
\end{aligned} − f ( s + 1 ) = n − s n + s − n 2 s + 1 f ( s ) + n − s s 2 s + 1 + n − s 1 x = 0 ∑ s − 1 ( x s ) ( ( n − 2 x + s ) f ( x + 1 ) + ( n + 2 x − s ) f ( x ) )
大概是这个东西吧。可能。
O ( n 2 ) O(n^2) O ( n 2 ) 。
我们设
F ( x ) = ∑ c P ( X = c ) x c F(x)=\sum_c P(X=c)x^c F ( x ) = ∑ c P ( X = c ) x c ,则
F ′ ( 1 ) = ∑ c P ( X = c ) ⋅ c = E ( X ) F'(1)=\sum_c P(X=c)\cdot c=\mathbb{E}(X) F ′ ( 1 ) = ∑ c P ( X = c ) ⋅ c = E ( X ) 。
先规定一些记号:
f i x ( i ) = ( 1 n ) i fix(i)=\left(\frac1n\right)^i f i x ( i ) = ( n 1 ) i ,
S ( i ) S(i) S ( i ) 表示一个长为
i i i 的字符串,
[ b o r d e r ( i ) ] [border(i)] [ b or d er ( i )] 指示
A [ 1 : i ] A[1:i] A [ 1 : i ] 是否是
A A A 的 border。
f i f_i f i 表示字符串延长到长度为
i i i 时结束的概率,
g i g_i g i 表示字符串长度为
i i i 时仍未结束的概率。
F ( x ) = ∑ i = 0 ∞ f i x i , G ( x ) = ∑ i = 0 ∞ g i x i F(x)=\sum_{i=0}^{\infty}f_ix^i,G(x)=\sum_{i=0}^{\infty}g_ix^i F ( x ) = ∑ i = 0 ∞ f i x i , G ( x ) = ∑ i = 0 ∞ g i x i 。
对于一个未结束的状态,增加一个字符要么结束要么未结束,即
g i = g i + 1 + f i + 1 g_i=g_{i+1}+f_{i+1} g i = g i + 1 + f i + 1 ,写成生成函数形式如下:
x G ( x ) + 1 = F ( x ) + G ( x ) xG(x)+1=F(x)+G(x) x G ( x ) + 1 = F ( x ) + G ( x )
然后整点小活。
( x G ( x ) ) ′ = F ′ ( x ) + G ′ ( x ) G ( x ) + x G ′ ( x ) = F ′ ( x ) + G ′ ( x ) F ′ ( 1 ) = G ( 1 ) (xG(x))'=F'(x)+G'(x)\\
G(x)+xG'(x)=F'(x)+G'(x)\\
F'(1)=G(1) ( x G ( x ) ) ′ = F ′ ( x ) + G ′ ( x ) G ( x ) + x G ′ ( x ) = F ′ ( x ) + G ′ ( x ) F ′ ( 1 ) = G ( 1 )
于是我们只需要求
G ( 1 ) G(1) G ( 1 ) 就好了。
对于一个
S ( i ) S(i) S ( i ) ,如果它还没结束,我们在后面拼上一个
A A A 便可以强制让它结束。保证拼上一个
A A A 需要提供
f i x ( m ) fix(m) f i x ( m ) 的概率。我们可以写出
f i x ( m ) ⋅ g i = f i + m fix(m)\cdot g_i=f_{i+m} f i x ( m ) ⋅ g i = f i + m 。
上面这个东西当然是错的。我们可能刚拼上
A A A 的一个前缀就凑出了一个
A A A 。这样的前缀也是
A A A 的后缀,即它是一个 border。
重写上面的式子:
f i x ( m ) g i = ∑ j = 1 m [ border ( j ) ] ⋅ f i x ( m − j ) f i + j fix(m)g_i=\sum_{j=1}^m[\operatorname{border}(j)]\cdot fix(m-j)f_{i+j} f i x ( m ) g i = ∑ j = 1 m [ border ( j )] ⋅ f i x ( m − j ) f i + j ,然后改成生成函数。
f i x ( m ) x m G ( x ) = ∑ i = 1 m [ border ( i ) ] ⋅ f i x ( m − i ) x m − i F ( x ) fix(m)x^mG(x)=\sum_{i=1}^m[\operatorname{border}(i)]\cdot fix(m-i)x^{m-i}F(x) f i x ( m ) x m G ( x ) = ∑ i = 1 m [ border ( i )] ⋅ f i x ( m − i ) x m − i F ( x )
f i x ( m ) G ( 1 ) = F ( 1 ) ∑ i = 1 m [ border ( i ) ] ⋅ f i x ( m − i ) G ( 1 ) = F ( 1 ) ∑ i = 1 m [ border ( i ) ] ⋅ n i fix(m)G(1)=F(1)\sum_{i=1}^m[\operatorname{border}(i)]\cdot fix(m-i)\\
G(1)=F(1)\sum_{i=1}^m[\operatorname{border}(i)]\cdot n^i f i x ( m ) G ( 1 ) = F ( 1 ) i = 1 ∑ m [ border ( i )] ⋅ f i x ( m − i ) G ( 1 ) = F ( 1 ) i = 1 ∑ m [ border ( i )] ⋅ n i
F ( 1 ) F(1) F ( 1 ) 即
∑ f i \sum f_i ∑ f i 显然为
1 1 1 。所以我们就做完了。
线性规划
感觉这个挺大件的,确实不是一天两天就能搞明白的东西(
一个线性规划问题的形式如下:给定
c 1 , c 2 , … , c n c_1,c_2,\dotsc,c_n c 1 , c 2 , … , c n ,你需要确定一组
x 1 , x 2 , … , x n x_1,x_2,\dotsc,x_n x 1 , x 2 , … , x n ,使其在满足另外给出的一些限制的前提下,
∑ j = 1 n c j x j \sum_{j=1}^nc_jx_j ∑ j = 1 n c j x j 最大(或最小)。
每组限制都是如下形式:
f i ( x 1 , x 2 , … , x n ) = ∑ j = 1 n a i j x j ≤ b j f_i(x_1,x_2,\dotsc,x_n)=\sum_{j=1}^na_{ij}x_j\leq b_{j} f i ( x 1 , x 2 , … , x n ) = ∑ j = 1 n a ij x j ≤ b j 。每组的
≤ \leq ≤ 也可以独立地换成
≥ \geq ≥ 或
= = = 。
标准形
任意一个线性规划问题都可以转化成下面的模型(标准形):
max ∑ j = 1 n c j x j s . t . ∑ j = 1 n a i j x j ≤ b i i ∈ [ 1 , m ] x j ≥ 0 j ∈ [ 1 , n ] \max \sum_{j=1}^nc_jx_j\\
\mathrm{s.t.}\quad\begin{aligned}
&\sum_{j=1}^n a_{ij}x_j\leq b_i\qquad&i\in[1,m]\\
&\ x_j\geq 0 &j\in[1,n]
\end{aligned} max j = 1 ∑ n c j x j s.t. j = 1 ∑ n a ij x j ≤ b i x j ≥ 0 i ∈ [ 1 , m ] j ∈ [ 1 , n ]
具体的转化方式是:
如果目标为求最小值,将所有 c j c_j c j 取反求最大。
如果限制是 ≥ \geq ≥ ,将两边同乘 − 1 -1 − 1 。
如果限制是 = = = 拆成 f i ≥ b i f_i\geq b_i f i ≥ b i 、− f i ≥ − b i -f_i\geq -b_i − f i ≥ − b i 即可。
如果某个 x j x_j x j 没有取值限制,将其替换成两个新变量的差 x p − x q x_p-x_q x p − x q 。
当然,我们也可以使用矩阵记号简写:
max c x s . t . A x ≤ b x ≥ 0 \max\quad\mathbf{cx}\\
\mathrm{s.t.}\quad\begin{aligned}
&\mathbf{A}\mathbf{x}\leq\mathbf{b}\\
&\mathbf{x}\geq\mathbf{0}
\end{aligned} max cx s.t. Ax ≤ b x ≥ 0
其中向量的不等号定义为偏序关系,即每一个分量都
≤ \leq ≤ (或
≥ \geq ≥ )。
对偶问题
考虑另一种刻画目标函数的方法:
我们不求解
x \mathbf{x} x ,而是新构造一组
y \mathbf{y} y ,使
y i y_i y i 作为约束
i i i ——
∑ j = 1 n a i j x j \sum_{j=1}^na_{ij}x_j ∑ j = 1 n a ij x j 的系数,试图用其刻画目标函数的极值。
此时,每个 y i y_i y i 各代表了一个约束。
相当于,原先
c j c_j c j 表示第
j j j 个产品要卖多少钱,需要消耗
a i j a_{ij} a ij 个原料
i i i 制造;现今直接将原料
i i i 定价
y i y_i y i 出售,寻求一个定价方案使得收益比全做成产品更值。
则我们需要有:
∀ x j , ∑ j = 1 n c j x j ≤ ∑ i = 1 m y i ∑ j = 1 n a i j x j \forall x_j,\quad\sum_{j=1}^nc_jx_j\leq \sum_{i=1}^m y_i\sum_{j=1}^na_{ij}x_j ∀ x j , j = 1 ∑ n c j x j ≤ i = 1 ∑ m y i j = 1 ∑ n a ij x j
即:
c j ≤ ∑ i = 1 m a i j y i j ∈ [ 1 , n ] c_j\leq \sum_{i=1}^ma_{ij}y_i\quad j\in[1,n] c j ≤ i = 1 ∑ m a ij y i j ∈ [ 1 , n ]
当
y i ≥ 0 y_i\geq 0 y i ≥ 0 时,目标函数的一个上界即为
∑ i = 1 m y i b i \sum_{i=1}^my_ib_i ∑ i = 1 m y i b i ,表示定价
y i y_i y i 时全卖掉最多赚这么多。
对于买家,你想在卖家不亏的前提下尽可能省钱地买下所有原料,故写了一个线性规划:
min ∑ i = 1 m b i y i s . t . ∑ i = 1 m a i j y i ≥ c j j ∈ [ 1 , n ] y i ≥ 0 i ∈ [ 1 , m ] \min\sum_{i=1}^m b_iy_i\\
\mathrm{s.t.}\quad\begin{aligned}
&\sum_{i=1}^ma_{ij}y_i\geq c_j\qquad &j\in[1,n]\\
&\ y_i\geq 0\qquad &i\in[1,m]
\end{aligned} min i = 1 ∑ m b i y i s.t. i = 1 ∑ m a ij y i ≥ c j y i ≥ 0 j ∈ [ 1 , n ] i ∈ [ 1 , m ]
也就是:
min b ⋅ y s . t . A ⊤ y ≥ c y ≥ 0 \min\quad\mathbf{b}\cdot\mathbf{y}\\
\mathrm{s.t.}\quad\begin{aligned}
&\mathbf{A}^\top\mathbf{y}\geq\mathbf{c}\\
&\mathbf{y}\geq0
\end{aligned} min b ⋅ y s.t. A ⊤ y ≥ c y ≥ 0
这两个问题互为对偶问题,对偶问题再对偶回去是原问题。
强对偶性
有定理保证,若原问题可行且有界,则对偶问题的最优解和原问题相等。
另有互松弛定理说明,设
x \bf{x} x 和
y \bf{y} y 为两问题的可行解,则其同时为最优解当且仅当:
∀ j ∈ [ 1 , n ] \forall j\in[1,n] ∀ j ∈ [ 1 , n ] ,有
x j = 0 x_j=0 x j = 0 或
∑ i = 1 m a i j y j = c j \sum_{i=1}^ma_{ij}y_j=c_j ∑ i = 1 m a ij y j = c j ;
∀ i ∈ [ 1 , m ] \forall i\in[1,m] ∀ i ∈ [ 1 , m ] ,有
y i = 0 y_i=0 y i = 0 或
∑ j = 1 n a i j x i = b i \sum_{j=1}^na_{ij}x_i=b_i ∑ j = 1 n a ij x i = b i 。
单纯形法
这部分我好像写的很烂,并且意义其实不大。
考虑一个东西叫做松弛型。 具体地,在每一个限制左边添加一个松弛变量
x i ′ ≥ 0 x'_i\geq 0 x i ′ ≥ 0 ,变为等式:
x i ′ = b i − ∑ j = 1 n a i j x j x'_i=b_i-\sum_{j=1}^{n}a_{ij}x_j x i ′ = b i − ∑ j = 1 n a ij x j
我们称等式左边的为基变量,包含在右边的为非基变量。显然基变量永远有
m m m 个,且和限制一一对应。初始时基变量就是所有的松弛变量。松弛变量和原有变量没有本质区别,下面
将 x i ′ x'_i x i ′ 重命名为 x n + i x_{n+i} x n + i 。
如果我们让所有非基变量
= 0 =0 = 0 ,则可以得到一个
基本解 :
x i ′ = b i x'_i=b_i x i ′ = b i 。单纯形法的思想是,通过替换基变量,不断优化基本解,尝试增大目标函数的取值,直到到达最优解。
(p.s. 如果
b i < 0 b_i<0 b i < 0 ,基本解其实是不合法的。但是没关系,后面可以通过简单的方法让所有
b i ≥ 0 b_i\geq 0 b i ≥ 0 ,这里不妨先假定如此。)
先考虑如何判断最优解。我们在算法进行过程中维护
c \mathbf{c} c 使得对于所有基变量
x j x_j x j ,
c j = 0 c_j=0 c j = 0 。初始时显然满足。如果对于所有非基变量
x j x_j x j ,
c j ≤ 0 c_j\leq 0 c j ≤ 0 ,则不管增加任何一个非基变量都不能让目标函数取值更大,所以此时的基本解即为最优解。
如果对于某个非基变量
x j x_j x j ,
c j > 0 c_j>0 c j > 0 ,则增加
x j x_j x j 可以使目标函数取值更大,考虑让它替代某个基变量。如果选择使用限制
i i i 让
x j x_j x j 成为基变量,则需要将限制
i i i 的所有系数(包括常数
b i b_i b i )乘上某个数
t t t ,然后再加减消元灭掉其它限制(
和目标函数 c \mathbf{c} c )中
x j x_j x j 的系数。其实这就是做了一套
初等行变换 ,和高消是一样的。注意这样操作以后
c \mathbf{c} c 可能会带常数,别不管。
“使用限制
i i i 让
x j x_j x j 成为基变量” 这个操作,下称
pivot ( i , j ) \operatorname{pivot}(i,j) pivot ( i , j ) 。
如果我们选择某个限制
i i i 去
pivot ( i , j ) \operatorname{pivot}(i,j) pivot ( i , j ) ,则
a i j a_{ij} a ij 必须
> 0 >0 > 0 。如果
= 0 =0 = 0 没法将系数归一,而
< 0 <0 < 0 会让
b i < 0 b_i<0 b i < 0 ,非可行解。如果可行,
pivot ( i , j ) \operatorname{pivot}(i,j) pivot ( i , j ) 后
x j x_j x j 会取到
b i a i j \frac{b_i}{a_{ij}} a ij b i 。在众多可行的
i i i 中,我们选取
b i a i j \frac{b_i}{a_{ij}} a ij b i 最小的限制进行
pivot \operatorname{pivot} pivot ,否则操作后比其更小的限制就不满足了。如果找不到可行限制则
x j x_j x j 可以无限增大,所以结果无界。
现在可以处理初始
b i < 0 b_i<0 b i < 0 的问题了。随机地找到一个
< 0 <0 < 0 的
b i b_i b i ,再随机地找到一个
< 0 <0 < 0 的
a i j a_{ij} a ij ,
pivot ( i , j ) \operatorname{pivot}(i,j) pivot ( i , j ) 即可。注意必须带有随机性,如果不随机可能会死循环,鬼知道为什么。这个方法非常的玄学。
由于神秘的原因,随便选择非基变量
pivot \operatorname{pivot} pivot 可能也会死循环。具体实现上,我们每次选取第一个
c j > 0 c_j>0 c j > 0 的
x j x_j x j 进行
pivot \operatorname{pivot} pivot 。
上面的东西很乱,梳理一下需要干什么:
初始化可行解:用
pivot \operatorname{pivot} pivot 处理掉
b i < 0 b_i<0 b i < 0 。
找到第一个
c j < 0 c_j<0 c j < 0 的
x j x_j x j 。如果找不到就说明当前是最优解,
c \mathbf{c} c 带的常数取相反数就是答案。
找到限制最紧的限制
i i i ,
pivot ( i , j ) \operatorname{pivot}(i,j) pivot ( i , j ) 。找不到则答案无界。
手模一组样例。
max x 1 + 14 x 2 + 6 x 3 s . t . { x 1 + x 2 + x 3 + x 4 = 4 x 1 + x 5 = 2 x 3 + x 6 = 3 3 x 2 + x 3 + x 7 = 6 x 1 , x 2 , x 3 , x 4 , x 5 , x 6 , x 7 ≥ 0 \max\quad x_1+14x_2+6x_3\\
\mathrm{s.t.}\quad\begin{cases}
x_1+x_2+x_3+x_4=4\\
x_1+x_5=2\\
x_3+x_6=3\\
3x_2+x_3+x_7=6\\
x_1,x_2,x_3,x_4,x_5,x_6,x_7\geq 0
\end{cases} max x 1 + 14 x 2 + 6 x 3 s.t. ⎩ ⎨ ⎧ x 1 + x 2 + x 3 + x 4 = 4 x 1 + x 5 = 2 x 3 + x 6 = 3 3 x 2 + x 3 + x 7 = 6 x 1 , x 2 , x 3 , x 4 , x 5 , x 6 , x 7 ≥ 0
初始:
x 1 x_1 x 1 x 2 x_2 x 2 x 3 x_3 x 3 x 4 x_4 x 4 x 5 x_5 x 5 x 6 x_6 x 6 x 7 x_7 x 7 b b b c 1 = 1 c_1=1 c 1 = 1 c 2 = 14 c_2=14 c 2 = 14 c 3 = 6 c_3=6 c 3 = 6 c 4 c_4 c 4 c 5 c_5 c 5 c 6 c_6 c 6 c 7 c_7 c 7 − a n s -ans − an s 1 1 1 1 1 1 1 1 1 1 1 1 4 4 4 1 1 1 1 1 1 2 2 2 1 1 1 1 1 1 3 3 3 3 3 3 1 1 1 1 1 1 6 6 6
将
x 1 x_1 x 1 变为基变量。限制
1 1 1 :
x 1 = 4 x_1=4 x 1 = 4 。限制
2 2 2 :
x 1 = 2 x_1=2 x 1 = 2 。
pivot ( 2 , 1 ) \operatorname{pivot}(2,1) pivot ( 2 , 1 ) 。
x 1 x_1 x 1 x 2 x_2 x 2 x 3 x_3 x 3 x 4 x_4 x 4 x 5 x_5 x 5 x 6 x_6 x 6 x 7 x_7 x 7 b b b c 1 c_1 c 1 c 2 = 14 c_2=14 c 2 = 14 c 3 = 6 c_3=6 c 3 = 6 c 4 c_4 c 4 c 5 = − 1 c_5=-1 c 5 = − 1 c 6 c_6 c 6 c 7 c_7 c 7 − a n s = − 2 -ans=-2 − an s = − 2 1 1 1 1 1 1 1 1 1 − 1 -1 − 1 2 2 2 1 1 1 1 1 1 2 2 2 1 1 1 1 1 1 3 3 3 3 3 3 1 1 1 1 1 1 6 6 6
将
x 2 x_2 x 2 变为基变量。限制
1 1 1 :
x 2 = 2 x_2=2 x 2 = 2 。限制
4 4 4 :
x 2 = 2 x_2=2 x 2 = 2 。
pivot ( 1 , 2 ) \operatorname{pivot}(1,2) pivot ( 1 , 2 ) 。
x 1 x_1 x 1 x 2 x_2 x 2 x 3 x_3 x 3 x 4 x_4 x 4 x 5 x_5 x 5 x 6 x_6 x 6 x 7 x_7 x 7 b b b c 1 c_1 c 1 c 2 c_2 c 2 c 3 = − 8 c_3=-8 c 3 = − 8 c 4 = − 14 c_4=-14 c 4 = − 14 c 5 = 13 c_5=13 c 5 = 13 c 6 c_6 c 6 c 7 c_7 c 7 − a n s = − 30 -ans=-30 − an s = − 30 1 1 1 1 1 1 1 1 1 − 1 -1 − 1 2 2 2 1 1 1 1 1 1 2 2 2 1 1 1 1 1 1 3 3 3 − 2 -2 − 2 − 3 -3 − 3 3 3 3 1 1 1 0 0 0
将
x 5 x_5 x 5 变为基变量。限制
2 2 2 :
x 5 = 2 x_5=2 x 5 = 2 。限制
4 4 4 :
x 5 = 0 x_5=0 x 5 = 0 。
pivot ( 4 , 5 ) \operatorname{pivot}(4,5) pivot ( 4 , 5 ) 。
x 1 x_1 x 1 x 2 x_2 x 2 x 3 x_3 x 3 x 4 x_4 x 4 x 5 x_5 x 5 x 6 x_6 x 6 x 7 x_7 x 7 b b b c 1 c_1 c 1 c 2 c_2 c 2 c 3 = 2 / 3 c_3=2/3 c 3 = 2/3 c 4 = − 1 c_4=-1 c 4 = − 1 c 5 c_5 c 5 c 6 c_6 c 6 c 7 = − 13 / 3 c_7=-13/3 c 7 = − 13/3 − a n s = − 30 -ans=-30 − an s = − 30 1 1 1 1 / 3 1/3 1/3 1 / 3 1/3 1/3 2 2 2 1 1 1 2 / 3 2/3 2/3 1 1 1 − 1 / 3 -1/3 − 1/3 2 2 2 1 1 1 1 1 1 3 3 3 − 2 / 3 -2/3 − 2/3 − 1 -1 − 1 1 1 1 1 / 3 1/3 1/3 0 0 0
将
x 3 x_3 x 3 变为基变量。限制
1 1 1 :
x 3 = 6 x_3=6 x 3 = 6 。限制
2 2 2 :
x 3 = 3 x_3=3 x 3 = 3 。限制
3 3 3 :
x 3 = 3 x_3=3 x 3 = 3 。
pivot ( 3 , 3 ) \operatorname{pivot}(3,3) pivot ( 3 , 3 ) 。
x 1 x_1 x 1 x 2 x_2 x 2 x 3 x_3 x 3 x 4 x_4 x 4 x 5 x_5 x 5 x 6 x_6 x 6 x 7 x_7 x 7 b b b c 1 = − 1 c_1=-1 c 1 = − 1 c 2 c_2 c 2 c 3 c_3 c 3 c 4 = − 2 c_4=-2 c 4 = − 2 c 5 c_5 c 5 c 6 c_6 c 6 c 7 = − 4 c_7=-4 c 7 = − 4 − a n s = − 32 -ans=-32 − an s = − 32 − 1 / 2 -1/2 − 1/2 1 1 1 − 1 / 2 -1/2 − 1/2 1 / 2 1/2 1/2 1 1 1 3 / 2 3/2 3/2 1 1 1 3 / 2 3/2 3/2 − 1 / 2 -1/2 − 1/2 3 3 3 − 3 / 2 -3/2 − 3/2 − 3 / 2 -3/2 − 3/2 1 1 1 1 / 2 1/2 1/2 0 0 0 1 1 1 1 1 1 2 2 2
done.下班。
exgcd
求解不定方程
a x + b y = c ax+by=c a x + b y = c 。首先如果
c ∤ ( a , b ) c\nmid (a,b) c ∤ ( a , b ) 无解,否则可以转为
a x + b y = ( a , b ) ax+by=(a,b) a x + b y = ( a , b ) 的一般问题上。
设
a x + b y = ( a , b ) ax+by=(a,b) a x + b y = ( a , b ) ,
b x ′ + ( a m o d b ) y ′ = ( b , a m o d b ) bx'+(a\bmod b)y'=(b,a\bmod b) b x ′ + ( a mod b ) y ′ = ( b , a mod b ) 。则因为
( a , b ) = ( b , a m o d b ) (a,b)=(b,a\bmod b) ( a , b ) = ( b , a mod b ) ,有:
a x + b y = b x ′ + ( a m o d b ) y ′ a x + b y = b x ′ + ( a − b ⌊ a b ⌋ ) y ′ a ( x − y ′ ) + b ( y − x ′ + ⌊ a b ⌋ y ′ ) = 0 \begin{aligned}
ax+by&=bx'+(a\bmod b)y'\\
ax+by&=bx'+(a-b\lfloor\frac{a}{b}\rfloor)y'\\
a(x-y')+b(y-x'+\lfloor\frac{a}{b}\rfloor y')&=0
\end{aligned} a x + b y a x + b y a ( x − y ′ ) + b ( y − x ′ + ⌊ b a ⌋ y ′ ) = b x ′ + ( a mod b ) y ′ = b x ′ + ( a − b ⌊ b a ⌋) y ′ = 0
所以
{ x = y ′ y = x ′ − ⌊ a b ⌋ y ′ \begin{cases}x=y'\\y=x'-\lfloor\frac{a}{b}\rfloor y'\end{cases} { x = y ′ y = x ′ − ⌊ b a ⌋ y ′ 。终止条件是
b = 0 b=0 b = 0 时有
{ x = 1 y = 0 \begin{cases}x=1\\y=0\end{cases} { x = 1 y = 0 。
exCRT
求解
{ x ≡ a 1 ( m o d m 1 ) x ≡ a 2 ( m o d m 2 ) ⋮ x ≡ a n ( m o d m n ) \begin{cases}x\equiv a_1\pmod {m_1}\\x\equiv a_2\pmod {m_2}\\\quad\vdots\\ x\equiv a_n\pmod {m_n}\end{cases} ⎩ ⎨ ⎧ x ≡ a 1 ( mod m 1 ) x ≡ a 2 ( mod m 2 ) ⋮ x ≡ a n ( mod m n ) 。
考虑一次将两个同余方程合并。
{ x ≡ a 1 ( m o d m 1 ) x ≡ a 2 ( m o d m 2 ) ⇔ x ≡ a ′ ( m o d lcm ( m 1 , m 2 ) ) \begin{cases}x\equiv a_1\pmod {m_1}\\x\equiv a_2\pmod {m_2}\end{cases}\Leftrightarrow x\equiv a'\pmod{\operatorname{lcm}(m_1,m_2)} { x ≡ a 1 ( mod m 1 ) x ≡ a 2 ( mod m 2 ) ⇔ x ≡ a ′ ( mod lcm ( m 1 , m 2 ))
其中
a ′ = a 1 + k 1 ⋅ m 1 = a 2 + k 2 ⋅ m 2 a'=a_1+k_1\cdot m_1=a_2+k_2\cdot m_2 a ′ = a 1 + k 1 ⋅ m 1 = a 2 + k 2 ⋅ m 2 ,这是 exgcd 需要做的事情。
Lucas
( n m ) m o d p = ( ⌊ n / p ⌋ ⌊ m / p ⌋ ) ⋅ ( n m o d p m m o d p ) m o d p \binom{n}{m}\bmod p=\binom{\lfloor n/p\rfloor}{\lfloor m/p\rfloor}\cdot \binom{n\bmod p}{m\bmod p}\bmod p ( m n ) mod p = ( ⌊ m / p ⌋ ⌊ n / p ⌋ ) ⋅ ( m mod p n mod p ) mod p
不要管它怎么来的。
FWT
设
U = [ 0 , n ) ∩ N U=[0,n)\cap \mathbb{N} U = [ 0 , n ) ∩ N ,给定两个集合幂级数
A , B : 2 U → R A,B:2^U\to R A , B : 2 U → R ,求
C S = ∑ T 1 ⊙ T 2 = S A T 1 B T 2 C_S=\sum_{T_1\odot T_2=S}A_{T_1}B_{T_2} C S = ∑ T 1 ⊙ T 2 = S A T 1 B T 2 。
⊙ \odot ⊙ 是某种位运算。暴力是
O ( 3 n ) O(3^n) O ( 3 n ) 的,很坏。
下面可能会用
N \mathbb{N} N 表示其二进制集合。
或
记数列
F ( X ) F(X) F ( X ) 为幂级数
X X X 的变换,满足
F ( X ) S = ∑ T ⊆ S X T F(X)_S=\sum_{T\sube S} X_T F ( X ) S = ∑ T ⊆ S X T 。可知,
F ( A ) S ⋅ F ( B ) S = F ( C ) S F(A)_S\cdot F(B)_S=F(C)_S F ( A ) S ⋅ F ( B ) S = F ( C ) S 。下面浅证一下。
F ( A ) S ⋅ F ( B ) S = ∑ P ⊆ S ∑ Q ⊆ S A P B Q = ∑ ( P ∪ Q ) ⊆ S A P B Q = F ( C ) S \begin{aligned}
F(A)_S\cdot F(B)_S&=\sum_{P\sube S}\sum_{Q\sube S}A_PB_Q\\
&=\sum_{(P\cup Q)\sube S}A_PB_Q\\
&=F(C)_S
\end{aligned} F ( A ) S ⋅ F ( B ) S = P ⊆ S ∑ Q ⊆ S ∑ A P B Q = ( P ∪ Q ) ⊆ S ∑ A P B Q = F ( C ) S
所以可以认为
F ( X ) F(X) F ( X ) 就是 FFT 里面类似点值表示法的东西,可以直接
O ( 2 n ) O(2^n) O ( 2 n ) 乘。只需考虑如何快速计算变换及逆变换即可。
考虑分治。记
L L L 和
R R R 分别为
X X X 的前半部分和后半部分(
L i = X i , R i = X 2 n − 1 + i , L_i=X_i,R_i=X_{2^{n-1}+i}, L i = X i , R i = X 2 n − 1 + i , )。假设已经有了
F ( L ) F(L) F ( L ) 和
F ( R ) F(R) F ( R ) ,则对于
L L L 的部分(
i ∈ [ 0 , 2 n − 1 ) i\in[0,2^{n-1}) i ∈ [ 0 , 2 n − 1 ) ),显然有
F ( X ) i = F ( L ) i F(X)_i=F(L)_i F ( X ) i = F ( L ) i 。对于
R R R 的部分(
i ∈ [ 2 n − 1 , 2 n ) i\in[2^{n-1},2^n) i ∈ [ 2 n − 1 , 2 n ) ),仅有
( i − 2 n − 1 ) ⊆ i (i-2^{n-1})\sube i ( i − 2 n − 1 ) ⊆ i 没在
R R R 内部算上。所以
F ( X ) i = F ( L ) i − 2 n − 1 + F ( R ) i − 2 n − 1 F(X)_i=F(L)_{i-2^{n-1}}+F(R)_{i-2^{n-1}} F ( X ) i = F ( L ) i − 2 n − 1 + F ( R ) i − 2 n − 1 。
换一个更直观的表示方法就是:
F ( X ) = merge ( F ( L ) , F ( L ) + F ( R ) ) F(X)=\operatorname{merge}(F(L),F(L)+F(R)) F ( X ) = merge ( F ( L ) , F ( L ) + F ( R ))
+ + + 就是按位加,
merge \operatorname{merge} merge 类似于字符串拼接。
当
X X X 的长度为
1 1 1 时,显然
F ( X ) = X F(X)=X F ( X ) = X 。所以我们就可以
O ( n log n ) O(n\log n) O ( n log n ) 计算变换了。而逆变换将
merge \operatorname{merge} merge 中的
+ + + 替换成
− - − 即可。显然。代码如下,
f=0 表示是逆变换。
CPP void OR (int a[],bool f) {
int x=(f?1 :-1 );
for (int w=2 ;w<=N;w<<=1 ){
int h=w>>1 ;
for (int i=0 ;i<N;i+=w) for (int j=0 ;j<h;j++){
a[i+j+h]+=x*a[i+j]; a[i+j+h]=(a[i+j+h]+P)%P;
}
}
}
如果严格按照上面说的,逆变换时 w 应该降序。但是考虑变换并不关心 w 的顺序,因为所有二进制位是完全等价的 。可以参考 SOSdp。
与
可以认为与本质上就是取反后或,或完了再取反回来,或者你认为
F ( X ) S = ∑ T ⊇ S X T F(X)_S=\sum_{T\supe S} X_T F ( X ) S = ∑ T ⊇ S X T 。总之
F ( X ) = merge ( F ( L ) + F ( R ) , F ( R ) ) F(X)=\operatorname{merge}(F(L)+F(R),F(R)) F ( X ) = merge ( F ( L ) + F ( R ) , F ( R )) 。逆变换变号。
异或
子集超集什么的好像不太行了。考虑一下上面这些东西的本质。
记
c ( S , T ) c(S,T) c ( S , T ) 为
X T X_T X T 对
F ( X ) S F(X)_S F ( X ) S 的贡献系数,即
F ( X ) S = ∑ T c ( S , T ) X T F(X)_S=\sum_Tc(S,T)X_T F ( X ) S = ∑ T c ( S , T ) X T 。因为我们需要
F ( A ) S ⋅ F ( B ) S = F ( C ) S F(A)_S\cdot F(B)_S=F(C)_S F ( A ) S ⋅ F ( B ) S = F ( C ) S ,所以需要有
c ( S , T 1 ) ⋅ c ( S , T 2 ) = c ( S , T 1 ⊙ T 2 ) c(S,T_1)\cdot c(S,T_2)=c(S,T_1\odot T_2) c ( S , T 1 ) ⋅ c ( S , T 2 ) = c ( S , T 1 ⊙ T 2 ) 。
然后
c ( S , T ) c(S,T) c ( S , T ) 需要能拆到每一位上,即
c ( i , j ) = ∏ k c ( S ∩ { k } , T ∩ { k } ) c(i,j)=\prod_{k} c(S\cap \{k\},T\cap \{k\}) c ( i , j ) = ∏ k c ( S ∩ { k } , T ∩ { k }) 。我们只需要给出一个矩阵
[ c ( 0 , 0 ) c ( 0 , 1 ) c ( 1 , 0 ) c ( 1 , 1 ) ] \begin{bmatrix}c(0,0)&c(0,1)\\c(1,0)&c(1,1)\end{bmatrix} [ c ( 0 , 0 ) c ( 1 , 0 ) c ( 0 , 1 ) c ( 1 , 1 ) ] 即可完成变换的构造。对于或卷积,矩阵即为
[ 1 0 1 1 ] \begin{bmatrix}1&0\\1&1\end{bmatrix} [ 1 1 0 1 ] 。然后根据什么魔法结论,逆变换的矩阵的元素同样可以拆位并且应该就是上面拆位后的变换矩阵的逆矩阵,对于或即为
[ 1 0 − 1 1 ] \begin{bmatrix}1&0\\-1&1\end{bmatrix} [ 1 − 1 0 1 ] 。
那我们可以构造出来异或的变换矩阵:
[ 1 1 1 − 1 ] \begin{bmatrix}1&1\\1&-1\end{bmatrix} [ 1 1 1 − 1 ] 。逆矩阵为
[ 1 2 1 2 1 2 − 1 2 ] \begin{bmatrix}\frac{1}{2}&\frac{1}{2}\\\frac{1}{2}&-\frac{1}{2}\end{bmatrix} [ 2 1 2 1 2 1 − 2 1 ] 。说人话就是:
F(X)=\operatorname{merge}(F(L)+F(R),F(L)-F(R))\\
X=\operatorname{merge}(\frac12(L+R),\frac12(L-R))\\
\end{cases}
另外,有时候也会用到一个实际上的意义:
F ( X ) i = ∑ x ( − 1 ) popcnt ( x ∩ i ) X x F(X)_i=\sum_{x}(-1)^{\operatorname{popcnt}(x\cap i)}X_x F ( X ) i = ∑ x ( − 1 ) popcnt ( x ∩ i ) X x
奇妙小性质?
可以认为
{ X } \{X\} { X } 构成一个交换幺环(是不是域?可能是,不知道)
( R 2 k , + , ⋅ ) (\mathbb{R}^{2^k},+,\cdot) ( R 2 k , + , ⋅ ) 。
A + B A+B A + B 定义为按位相加,
A ⋅ B A\cdot B A ⋅ B 定义为
⊙ \odot ⊙ 运算卷积。其加法单位元为
0 \mathbf{0} 0 ,乘法单位元为
[ 1 , 0 , 0 , ⋯ , 0 ] [1,0,0,\cdots,0] [ 1 , 0 , 0 , ⋯ , 0 ] 。
{ Y } \{Y\} { Y } 同样构成一个交换幺环
( R 2 k , + , ⋅ ) (\mathbb{R}^{2^k},+,\cdot) ( R 2 k , + , ⋅ ) ,不同的是
A ⋅ B A\cdot B A ⋅ B 定义为按位乘。这里的乘法单位元为
[ 1 , 1 , 1 , ⋯ , 1 ] [1,1,1,\cdots,1] [ 1 , 1 , 1 , ⋯ , 1 ] 。
F F F 是一个
{ X } → { Y } \{X\}\rightarrow \{Y\} { X } → { Y } 的线性映射。
F F F 的线性性是显然的,
F ( X ) F(X) F ( X ) 的每一位都是
X i X_i X i 的线性组合,而且我们上面甚至已经开始对
F F F 的转移矩阵分析性质了。
还注意到,
Y Y Y 在加法和乘法中,每一位是不相关的。(不知道该怎么准确地描述(
那这看起来很帅,可以整点活了。
记
C i = [ i ∈ A ] C_i=[i\in A] C i = [ i ∈ A ] ,则要计算的即为
∑ i = 1 n ∑ j ≠ 0 ( C i ) j \sum_{i=1}^n\sum_{j\neq 0}(C^i)_j ∑ i = 1 n ∑ j = 0 ( C i ) j 。暴力的做法是求
∑ i = 1 n F ( C i ) \sum_{i=1}^n F(C^i) ∑ i = 1 n F ( C i ) ,然后再
F − 1 F^{-1} F − 1 回去。由于线性性,这个东西就等于
∑ i = 1 n ( F ( C ) ) i \sum_{i=1}^n(F(C))^i ∑ i = 1 n ( F ( C ) ) i 。然后我们可以把
F ( C ) F(C) F ( C ) 每一位拆开来考虑。记
F ( C ) j = a j F(C)_j=a_j F ( C ) j = a j ,那就变成了对于每一个分量
j j j 求
∑ i = 1 n a j i \sum_{i=1}^na_j^i ∑ i = 1 n a j i 。那这个东西岂不是唐爆了。算完以后
F − 1 F^{-1} F − 1 回去就做完了。
记
M = 16 M=16 M = 16 ,总复杂度为
( M + log n ) 2 M (M+\log n)2^M ( M + log n ) 2 M 。
啊啊啊啊。这鬼东西真的让我啃了两天。长脑子了。
给定
n n n 个整数,每个数
∈ [ 0 , 2 m ) \in[0,2^m) ∈ [ 0 , 2 m ) 。对于每一个
c ∈ [ 0 , m ] c\in[0,m] c ∈ [ 0 , m ] ,求选出一个子集,使得其异或和的
popcnt \operatorname{popcnt} popcnt 恰好
= c =c = c 的方案数。
n ≤ 2 × 10 5 , m ≤ 53 n\leq 2\times 10^5,m\leq 53 n ≤ 2 × 1 0 5 , m ≤ 53 。
记这些数的线性基
A A A 张成的线性空间为
span ( A ) \operatorname{span}(A) span ( A ) ,
A A A 的大小为
rank ( A ) \operatorname{rank}(A) rank ( A ) 。
如果只选线性基里的数,那么每个
x ∈ span ( A ) x\in \operatorname{span}(A) x ∈ span ( A ) 被异或出来的方案数显然为
1 1 1 。但是现在多了
n − rank ( A ) n-\operatorname{rank}(A) n − rank ( A ) 个
0 0 0 可以选,所以对于任意的
x x x ,答案
× 2 n − rank ( A ) \times{2^{n-\operatorname{rank}(A)}} × 2 n − rank ( A ) 即可。
那我们只需要对于一个线性基考虑原问题就好了。
首先我们有一个不用脑子的做法。我们枚举整个线性空间可以得到答案。时间复杂度为
O ( 2 rank ( A ) ) O(2^{\operatorname{rank}(A)}) O ( 2 rank ( A ) ) 。下面再给出一种时间复杂度为
O ( 2 m − rank ( A ) ) O(2^{m-\operatorname{rank}(A)}) O ( 2 m − rank ( A ) ) 的做法,则可以通过选择两个解法中更快的那个做到
O ( 2 m / 2 ) O(2^{m/2}) O ( 2 m /2 ) 。
记集合幂级数
F ( A ) = ∑ x ∈ span ( A ) z x F(A)=\sum_{x\in\operatorname{span}(A)} z^x F ( A ) = ∑ x ∈ span ( A ) z x ,多项式
P ( A ) = ∑ x ∈ span ( A ) z popcnt ( x ) P(A)=\sum_{x\in\operatorname{span}(A)} z^{\operatorname{popcnt}(x)} P ( A ) = ∑ x ∈ span ( A ) z popcnt ( x ) 。以下指数
x x x 的加法均为
⊕ \oplus ⊕ ,集合幂级数的乘法直接遵循分配律(下面用
∗ * ∗ ),FWT 的乘法是按位乘(下面用
⋅ \cdot ⋅ )
对于一个幂级数
X X X ,下面不太区分
[ z i ] X [z^i]X [ z i ] X 和
X i X_i X i 两种表达,哪个好看用哪个。
我们想干的事情是,对于每一个
c ∈ [ 0 , m ] c\in[0,m] c ∈ [ 0 , m ] ,求
[ z c ] P ( A ) [z^c]P(A) [ z c ] P ( A ) 。
再记
G c = ∑ x z x [ popcnt ( x ) = c ] G_c=\sum_{x}z^x[\operatorname{popcnt}(x)=c] G c = ∑ x z x [ popcnt ( x ) = c ] ,则:
[ z c ] P ( A ) = [ z 0 ] F ( A ) ∗ G c = [ z 0 ] IFWT ( FWT ( F ( A ) ) ⋅ FWT ( G c ) ) [z^c]P(A)=[z^0]F(A)* G_c=[z^0]\operatorname{IFWT}(\operatorname{FWT}(F(A))\cdot \operatorname{FWT}(G_c)) [ z c ] P ( A ) = [ z 0 ] F ( A ) ∗ G c = [ z 0 ] IFWT ( FWT ( F ( A )) ⋅ FWT ( G c ))
这个式子我们现在还不会处理,别急,先做点其他的。屠龙要从第一刀开始切。
对于
x ∈ span ( A ) , z x ∗ F ( A ) = F ( A ) x\in\operatorname{span}(A),z^x* F(A)=F(A) x ∈ span ( A ) , z x ∗ F ( A ) = F ( A ) 。因为
∀ S ∈ span ( A ) , x + S ∈ span ( A ) \forall S\in \operatorname{span}(A),x+S\in\operatorname{span}(A) ∀ S ∈ span ( A ) , x + S ∈ span ( A ) 且显然
S S S 和
x + S x+S x + S 构成双射。
所以
F ( A ) ∗ F ( A ) = 2 rank ( A ) F ( A ) FWT ( F ( A ) ) ⋅ FWT ( F ( A ) ) = 2 rank ( A ) FWT ( F ( A ) ) F(A)*F(A)=2^{\operatorname{rank}(A)}F(A)\\
\operatorname{FWT}(F(A))\cdot \operatorname{FWT}(F(A))=2^{\operatorname{rank}(A)}\operatorname{FWT}(F(A)) F ( A ) ∗ F ( A ) = 2 rank ( A ) F ( A ) FWT ( F ( A )) ⋅ FWT ( F ( A )) = 2 rank ( A ) FWT ( F ( A ))
对每一位分别考虑。解方程
x 2 = 2 rank ( A ) x x^2=2^{\operatorname{rank}(A)}x x 2 = 2 rank ( A ) x ,得:
[ z i ] FWT ( F ( A ) ) ∈ { 0 , 2 rank ( A ) } [z^i]\operatorname{FWT}(F(A))\in\{0,2^{\operatorname{rank}(A)}\} [ z i ] FWT ( F ( A )) ∈ { 0 , 2 rank ( A ) }
发掘一下
FWT \operatorname{FWT} FWT 的性质。记
< i , j > = popcnt ( i & j ) \left<i,j\right>=\operatorname{popcnt}(i\&j) ⟨ i , j ⟩ = popcnt ( i & j ) ,则:
[ z i ] FWT ( X ) = ∑ x ( − 1 ) < i , x > X x [z^i]\operatorname{FWT}(X)=\sum_{x}(-1)^{\left<i,x\right>}X_x [ z i ] FWT ( X ) = ∑ x ( − 1 ) ⟨ i , x ⟩ X x
把
X X X 换成
F ( A ) F(A) F ( A ) 可得:
[ z i ] FWT ( F ( A ) ) = ∑ x ∈ span ( A ) ( − 1 ) < i , x > [z^i]\operatorname{FWT}(F(A))=\sum_{x\in\operatorname{span}(A)}(-1)^{\left<i,x\right>} [ z i ] FWT ( F ( A )) = ∑ x ∈ span ( A ) ( − 1 ) ⟨ i , x ⟩
发现好像当且仅当
∀ x , < i , x > ≡ 0 ( m o d 2 ) \forall x,\left<i,x\right>\equiv0\pmod 2 ∀ x , ⟨ i , x ⟩ ≡ 0 ( mod 2 ) 的时候
[ z i ] [z^i] [ z i ] 项才能取到
2 span ( A ) 2^{\operatorname{span}(A)} 2 span ( A ) ,否则根据上面的发现,只能是
0 0 0 。
我们还知道,
( i & x ) ⊕ ( j & x ) = ( i ⊕ j ) & x (i\& x)\oplus(j\& x)=(i\oplus j)\&x ( i & x ) ⊕ ( j & x ) = ( i ⊕ j ) & x 。而再推一步就一下发现
< i , x > + < j , x > \left<i,x\right>+\left<j,x\right> ⟨ i , x ⟩ + ⟨ j , x ⟩ 和
< i ⊕ j , x > \left<i\oplus j,x\right> ⟨ i ⊕ j , x ⟩ 的奇偶性是一样的。所以:
[ z i ] FWT ( F ( A ) ) = { 2 rank ( A ) ∀ x ∈ A , < i , x > ≡ 0 ( m o d 2 ) 0 otherwise. [z^i]\operatorname{FWT}(F(A))=\begin{cases}2^{\operatorname{rank}(A)}&\forall x\in A,\left<i,x\right>\equiv0\pmod 2\\0&\text{otherwise.}\end{cases} [ z i ] FWT ( F ( A )) = { 2 rank ( A ) 0 ∀ x ∈ A , ⟨ i , x ⟩ ≡ 0 ( mod 2 ) otherwise.
我们考虑找一个新的极大线性基
B B B ,使得
∀ x ∈ span ( A ) , y ∈ span ( B ) , < x , y > ≡ 0 ( m o d 2 ) \forall x\in\operatorname{span}(A),y\in\operatorname{span}(B),\left<x,y\right>\equiv 0\pmod 2 ∀ x ∈ span ( A ) , y ∈ span ( B ) , ⟨ x , y ⟩ ≡ 0 ( mod 2 ) 。可知
rank ( B ) = m − rank ( A ) \operatorname{rank}(B)=m-\operatorname{rank}(A) rank ( B ) = m − rank ( A ) 。
根据前面说的一车话,
2 rank ( A ) F ( B ) = FWT ( F ( A ) ) 2^{\operatorname{rank}(A)}F(B)=\operatorname{FWT}(F(A)) 2 rank ( A ) F ( B ) = FWT ( F ( A ))
怎么求这个
B B B ?考虑将
A A A 消元后,和
B B B 的基向量组成的矩阵。这个矩阵应该是关于主对角线对称的。反正能求就是了。
关键的部分差不多了,剩下一些杂碎。先处理
[ z i ] FWT ( G c ) [z^i]\operatorname{FWT}(G_c) [ z i ] FWT ( G c ) 。
展开:
[ z i ] FWT ( G c ) = ∑ popcnt ( x ) = c ( − 1 ) < i , x > = ∑ c n t = 0 m ( − 1 ) c n t ( popcnt ( i ) c n t ) ( m − popcnt ( i ) c − c n t ) [z^i]\operatorname{FWT}(G_c)=\sum_{\operatorname{popcnt}(x)=c}(-1)^{\left<i,x\right>}=\sum_{cnt=0}^m(-1)^{cnt}\binom{\operatorname{popcnt}(i)}{cnt}\binom{m-\operatorname{popcnt}(i)}{c-cnt} [ z i ] FWT ( G c ) = ∑ popcnt ( x ) = c ( − 1 ) ⟨ i , x ⟩ = ∑ c n t = 0 m ( − 1 ) c n t ( c n t popcnt ( i ) ) ( c − c n t m − popcnt ( i ) )
发现只和
popcnt ( i ) \operatorname{popcnt}(i) popcnt ( i ) 与
c c c 有关,可以
O ( m 3 ) O(m^3) O ( m 3 ) 预处理。记其为
g ( c , c n t ) g(c,cnt) g ( c , c n t ) 。
最后可以对外层的
IFWT \operatorname{IFWT} IFWT 下手了。根据计算
FWT \operatorname{FWT} FWT 的方法可以简单地得到:
IFWT ( X ) = 2 − m FWT ( X ) \operatorname{IFWT}(X)=2^{-m}\operatorname{FWT}(X) IFWT ( X ) = 2 − m FWT ( X )
然后找性质,可以得知:
[ z 0 ] FWT ( X ) = ∑ x X x [z^0]\operatorname{FWT}(X)=\sum_{x}X_x [ z 0 ] FWT ( X ) = ∑ x X x
我们爆搜出
P ( B ) P(B) P ( B ) ,则把上面一坨东西全部展开即有:
[ z c ] P ( A ) = 2 rank ( A ) − m ∑ d = 0 m [ z d ] P ( B ) g ( c , d ) [z^c]P(A)=2^{\operatorname{rank}(A)-m}\sum_{d=0}^m[z^d]P(B)g(c,d) [ z c ] P ( A ) = 2 rank ( A ) − m ∑ d = 0 m [ z d ] P ( B ) g ( c , d )
复杂度大头显然是
P ( B ) P(B) P ( B ) 的爆搜,即
2 rank ( B ) = 2 m − rank ( A ) 2^{\operatorname{rank}(B)}=2^{m-\operatorname{rank}(A)} 2 rank ( B ) = 2 m − rank ( A ) 。
所以这道题我们就杀了。
二项式反演
完蛋。这个东西我好像真不是很会。
g n = ∑ i = 0 n ( n i ) f i ⟺ f n = ∑ i = 0 n ( − 1 ) n − i ( n i ) g i g_n=\sum_{i=0}^n\binom{n}{i}f_i\Longleftrightarrow f_n=\sum_{i=0}^n(-1)^{n-i}\binom{n}{i}g_i g n = ∑ i = 0 n ( i n ) f i ⟺ f n = ∑ i = 0 n ( − 1 ) n − i ( i n ) g i
g k = ∑ i = k ( i k ) f i ⟺ f k = ∑ i = k ( − 1 ) i − k ( i k ) g i g_k=\sum_{i=k}\binom{i}{k} f_i\Longleftrightarrow f_k=\sum_{i=k}(-1)^{i-k}\binom{i}{k}g_i g k = ∑ i = k ( k i ) f i ⟺ f k = ∑ i = k ( − 1 ) i − k ( k i ) g i
第一个式子一般用作 “最多” 和 “恰好” 的转换。第二个式子一般用作 “至少” 和 “恰好” 的转换。
给定一个 DAG,对于每条边
( u , v ) (u,v) ( u , v ) 保证
u < v u<v u < v 。对于
∀ x ∈ [ 0 , k ] \forall x\in[0,k] ∀ x ∈ [ 0 , k ] ,求随机割掉
x x x 条边后,再随机选择原图上一条从
1 1 1 到
n n n 的路径,其未被割断的概率。
先考虑枚举
x x x 。将概率转化为方案数,即计算有多少种选择割掉的边和经过的路径的方案数,使得路径没有被割掉。计算出来后再除以总方案数
( m x ) ⋅ Path ( 1 , n ) \binom{m}{x}\cdot\operatorname{Path}(1,n) ( x m ) ⋅ Path ( 1 , n ) 即可。
对这个进行容斥。记
f i f_i f i 表示路径恰好经过
i i i 条割掉的边的方案数,
g i g_i g i 表示钦定路径上经过
i i i 条割掉的边的方案数。则:
g i = ∑ j = i ( j i ) f j g_i=\sum_{j=i}^{}\binom{j}{i}f_j g i = ∑ j = i ( i j ) f j
我们需要的是
f 0 = ∑ i = 0 ( − 1 ) i g i f_0=\sum_{i=0}(-1)^ig_i f 0 = ∑ i = 0 ( − 1 ) i g i 。那可以考虑怎么算一下
g i g_i g i 了。
记
h i h_i h i 表示在路径上随便钦定
i i i 条边的方案数,则
g i = ( m − i x − i ) h i g_i=\binom{m-i}{x-i}h_i g i = ( x − i m − i ) h i 。(已经钦定了的必须选,然后在剩下的边里随便选一些。)
因为
i ≥ x i\geq x i ≥ x 时系数
= 0 =0 = 0 ,且
x ≤ k x\leq k x ≤ k 所以我们只需要预处理出
h 0 , h 1 , ⋯ h k h_0,h_1,\cdots h_k h 0 , h 1 , ⋯ h k 就做完了。这是简单的。
其实算是一个trick。其实自己发明出来过。
考虑如何求一串多项式连乘
∏ i = 1 n F i \prod_{i=1}^{n}F_i ∏ i = 1 n F i 。你可以直接上
F F T \mathrm{FFT} FFT ,但是如果每个多项的次数都很小那就没什么优化。那么考虑分治求解,记
P ( 1 , n ) P(1,n) P ( 1 , n ) 为我们需要的,对于每个
P ( l , r ) P(l,r) P ( l , r ) 分解为
P ( l , m i d ) × P ( m i d + 1 , r ) P(l,mid)\times P(mid+1,r) P ( l , mi d ) × P ( mi d + 1 , r ) 求解即可。
O ( n log 2 n ) O(n\log^2 n) O ( n log 2 n ) 。
将答案写成多项式的形式,则:
F k = x + ∏ v ∈ son ( k ) F v F_k=x+\prod_{v\in \text{son}(k)}F_v F k = x + ∏ v ∈ son ( k ) F v
直接算是
O ( n 2 log n ) O(n^2\log n) O ( n 2 log n ) ,好像还完全不够。
因为时间开销和子树大小严格相关,考虑树剖,进行类似于动态 dp 的合并轻儿子贡献。记节点
k k k 的重儿子为
h ( k ) h(k) h ( k ) ,
G k = ∏ v ∈ son ( k ) ∧ v ≠ h ( k ) F v G_k=\prod_{v\in\operatorname{son}(k)\land v\neq h(k)}F_v G k = ∏ v ∈ son ( k ) ∧ v = h ( k ) F v ,则改写成:
F k = x + F h ( k ) G k F_k=x+F_{h(k)}G_k F k = x + F h ( k ) G k
然后考虑将后面的
F h ( k ) F_{h(k)} F h ( k ) 展开,本质上是将一整条重链拉出来考虑。
F k = x + ( x + ( x + ( ⋯ ) G h 2 ( k ) ) G h ( k ) ) G k F_k=x+(x+(x+(\cdots)G_{h^2(k)})G_{h(k)})G_k F k = x + ( x + ( x + ( ⋯ ) G h 2 ( k ) ) G h ( k ) ) G k
类似于秦九韶的形式。分别将
G k , G h ( k ) , G h 2 ( k ) , ⋯ G h s − 1 ( k ) G_k,G_{h(k)},G_{h^2(k)},\cdots G_{h^{s-1}(k)} G k , G h ( k ) , G h 2 ( k ) , ⋯ G h s − 1 ( k ) 重新标号为
G 1 , G 2 , G 3 , ⋯ G s G_{1},G_{2},G_{3},\cdots G_{s} G 1 , G 2 , G 3 , ⋯ G s ,把原式完全展开。
F k = x + ∏ i = 1 s G i + ∑ i = 1 s − 1 x ∏ j = 1 i G j F_k=x+\prod_{i=1}^{s}G_i+\sum_{i=1}^{s-1}x\prod_{j=1}^{i}G_j F k = x + ∏ i = 1 s G i + ∑ i = 1 s − 1 x ∏ j = 1 i G j
这个怎么做。中间的那项可以分治 FFT,后面一项其实也可以。记
A ( l , r ) = ∏ i = l r G i A(l,r)=\prod_{i=l}^{r}G_i A ( l , r ) = ∏ i = l r G i ,
B ( l , r ) = ∑ i = 1 r − 1 x ∏ j = l i G j B(l,r)=\sum_{i=1}^{r-1}x\prod_{j=l}^{i}G_j B ( l , r ) = ∑ i = 1 r − 1 x ∏ j = l i G j ,则:
A(l,r)&=A(l,m)A(m+1,r)\\
B(l,r)&=B(l,m)+B(m+1,r)A(l,m)+A(l,m)
\end{aligned}$$
那这样你就可以 $O(sz\log^2 n)$ 地做后面的东西了。
然后因为重剖了,所以所有链顶的 $sz$ 之和为 $O(n\log n)$。共 $O(n\log^3 n)$。
## [计数](https://www.luogu.com.cn/problem/AT_agc065_f)
写得太破防了,不太想放到做题记录里。
***
有人已经证了,这种图仅会由以下方式生成:
1. 在图上建立若干个偶环/两点一线,
2. 选择一些连通块,每个连通块里抽出来一个点,连成一个点双
3. 重复 2. 直到全图连通。
记 $g_n$ 为 $n$ 个点的带标号无向连通图个数,则根据 $1$ 所在的连通块讨论可得:
$$g_n=2^{\frac{1}{2}n(n-1)}-\sum_{i=1}^{n-1}\binom{n-1}{i-1}2^{\frac{1}{2}(n-i)(n-i-1)}g_i$$
然后记 $f_{n,m}$ 表示 $n$ 个点的连通图恰有 $m$ 个点双的方案数。$m\neq 1$ 时考虑枚举点双的大小拼起来:
$$f_{n,m}=\frac{1}{m!}\sum_{a_1+a_2+\cdots+a_m=n-1}\binom{n-1}{a_1,a_2,\cdots,a_m}n^{m-1}\prod_{i=1}^{m}f_{a_i+1,1}$$
考虑建出圆方树,则除了 $\Large\bigcirc\normalsize\kern{-9.5px}1\kern{5.5px}$ 以外,每一个方点的儿子数都是 $a_i-1$。将若干个连通块连起来的方案数根据 $\mathrm{Pr\ddot{u}fer}$ 序列的结论是 $n^{m-2}\prod a_i$,但是因为我们实际上使用方点而不是真实的圆点连边,所以方案数正好把连乘除下去了。然后需要考虑重新进行编号分配,就是那个组合数。$1/m!$ 是消除我们对点双钦定的顺序。
然后 $f_{n,1}=g_n-\sum_{i=2}^{n}f_{n,i}$。
考虑算答案。枚举放了 $k$ 个结构,第 $i$ 个结构的大小为 $b_i\ \ (\sum b_i=n)$,用 $m$ 个 $|S_i|\equiv 0\pmod2$ 的 DCC 将其粘起来。则有:
$$ans=\color{#53B396}{\binom{n}{b_1,b_2,\cdots,b_k}\frac{1}{k!}\prod_{i=1}^{k}\operatorname{calc}(b_i)\cdot b_i}\cdot \color{#9148AB}{\frac{1}{m!}\sum_{a_1+a_2+\cdots+a_m=k-1}n^{m-1}\binom{k-1}{a_1,a_2,\cdots,a_m}\prod_{i=1}^{m}f_{a_i+1,1}}$$
青色的部分在给所有结构标号。首先给所有环分配统一标号,然后使用 $\operatorname{calc}(x)=\begin{cases}1&x=2\\(x-1)!/2&\text{otherwise.}\end{cases}$ 计算一个环内部的不同标号方式。然后因为后面紫色部分不关心环的顺序,所以 $\div\ k!$ 消去环本身的编号。$b_i$ 为每个环钦定了一个 “根”,作用是规定这个环用哪个点挂到其在圆方树的父亲上。即使是 $\Large\bigcirc\normalsize\kern{-9.5px}1\kern{5.5px}$ 也需要一个,这是因为 $\mathrm{Pr\ddot{u}fer}$ 结论后面带一个 $\prod b_i$,只是因为其它连通块钦定方点连边抵消了。$\Large\bigcirc\normalsize\kern{-9.5px}1\kern{5.5px}$ 是自己连边不能抵消。
紫色部分在计算粘连方式,大体和上面一样,只不过我们现在需要具体考虑每个方点向哪个环上的点连边,所以底数是 $n$。
可以整理一下:
$$ans=\frac{n!}{k!}\frac{n^{m-1}}{k^{m-1}}f_{k,m}\prod_{i=1}^{k}\frac{1+[b_i=2]}{2}$$
虽然式子其实没问题但是应该没人处理了 $f_{1,0}$,所以整张图就是环的情况最好特判。
## 牛顿迭代法
给定二元函数 $G(x,y)$,求 $f(x)$ 使得 $G(x,f(x))\equiv 0\pmod {x^n}$。
假设已经获得了 $f(x)$ 的前 $k$ 项 $f_k(x)$,考虑怎么获得前 $2k$ 项。
将 $G(x,y)$ 对 $y$ 在 $y=f_k(x)$ 处泰勒展开:
$$G\big(x,f(x)\big)=\sum_{i=0}^{\infty}\frac{\partial^i}{\partial y^i}G\big(x,f_k(x)\big)\cdot \frac{\big(f(x)-f_k(x)\big)^i}{i!}\equiv 0\pmod{x^{2k}}$$
$\big(f(x)-f_k(x)\big)$ 的最低非 $0$ 项次数 $\geq k$,故 $\forall i\geq 2,\big(f(x)-f_k(x)\big)^i\equiv0\pmod {x^{2k}}$。故:
$$G\big(x,f_k(x)\big)+\frac{\partial G}{\partial y}\big(x,f_k(x)\big)\cdot \big(f(x)-f_k(x)\big)\equiv0\pmod {x^{2k}}$$
变形可得:
$$f_{2k}(x)\equiv f_k(x)-\frac{G\big(x,f_k(x)\big)}{\frac{\partial G}{\partial y}\big(x,f_k(x)\big)}\pmod {x^{2k}}$$
于是只要手动获得一下 $f_1(x)$ 就可以递推了。时间复杂度一般是 $T(n)=T(n/2)+O(n\log n)=O(n\log n)$。
### 多项式求逆
设 $G(x,y)=\dfrac{1}{y}-h(x)$,化简得:
$$f_{2k}(x)\equiv 2f_k(x)-f_k^2(x)h(x)\pmod {x^{2k}}$$
## 多项式 exp
设 $G(x,y)=\ln y-h(x)$,化简得:
$$f_{2k}(x)\equiv f_k(x)(1-\ln f_k(x)+h(x))\pmod {x^{2k}}$$
$\ln f(x)=\int \frac{f'(x)}{f(x)}\mathrm{d}x$,每次用多项式求逆现算,还是 $O(n\log n)$。注意这里 $\ln f_k(x)$ 要算到 $2k$ 项。