标签:algorithm math mp3 += deb algo list printf 转换
目录
参考博客-1
zjt 的博客地址挂掉了,所以这里暂时没有。
我们知道,有时候可以使用牛顿迭代法求 G(F(x)) = 0 的某个根 F0(x)。
但是有时候函数 G 是一个非常复杂的多项式。
这个时候使用牛顿迭代法,每一层的复杂度可能要做 n 次乘法与加法(n 为当前多项式长度)。
于是时间复杂度将会急剧退化,甚至不如暴力。
但假如 G(F(x)) = x,其中 G, F 都为多项式函数,在已知 F 的情况下,我们可以使用拉格朗日反演求出 G 的某一项。
同时,假如 G(F(x)) = H(x),我们也可以使用扩展拉格朗日反演求出在已知 F 的情况下求 G 的某一项。
形式化地给出问题:
若已知多项式函数 F(x),且复合函数 G(F(x)) 满足 G(F(x)) = x,求 G(x) 的第 n 项。
此时可以看作是 x 在复合函数 G·F 的作用下保持不变,所以我们称 G 为 F 的复合逆(复合运算下的逆元)。
不过对于 F,G 是有要求的:要求两个多项式常数项为 0 且一次项系数不为 0。不然复合逆是不存在的。
有一个结论:若 F 为 G 的复合逆,则 G 为 F 的复合逆。即 F(G(x)) = x 意味着 G(F(x)) = x。
所以:我们称 F 与 G 互为复合逆。
这个技巧可以在已知 F 的情况下,将 F(G(x)) = x 转为 G(F(x)) = x,有助于我们解决一些题目。
开始描述解法:
令 \(G(x) = \sum_{i=1}a_i*x^i\),那么就有:
\[\sum_{i=1}a_i*F^i(x) = x\]
两边求导,得到:
\[\sum_{i=1}i*a_i*F^{i-1}(x)*F'(x) = 1\]
为了得到 G 的第 n 项 \(a_n\),我们将两边除以 \(F^n(x)\),得到:
\[\sum_{i=1}i*a_i*F^{i-n-1}(x)*F'(x) = \frac{1}{F^n(x)}\]
(这个地方注意一点:我们此时函数并不一定为多项式——我们将函数推广到了负整数幂)
对于幂函数 \(F^a(x)\),当 a ≠ 0 时,它的导数为 \(a*F^{a-1}(x)*F'(x)\)。
反过来,当 b ≠ -1 时,\(F^{b}(x)*F'(x)\) 的积分为 \(\frac{F^{b+1}(x)}{b+1}\)。
于是对上面的式子进行进一步转换:
\[n*a_n*\frac{F'(x)}{F(x)}+\sum_{i=1}^{i\not = n}\frac{i}{i-n}*a_i*(F^{i-n}(x))' = \frac{1}{F^n(x)}\]
等式左右对应项系数相等。假如我们取 \(x^{-1}\) 的系数,则等式左边只留下 \(n*a_n*\frac{F'(x)}{F(x)}\) 这一部分(另一部分因为幂函数求导不会出现 \(x^{x-1}\) 这一项)。
于是就有:
\[[x^{-1}](n*a_n*\frac{F'(x)}{F(x)}) = [x^{-1}](\frac{1}{F^n(x)})\]
其中\([x^{-1}]\) 表示取 \(x^{-1}\) 的系数。
然后我们考虑一下 \(\frac{F'(x)}{F(x)}\) 中 \(x^{-1}\) 的系数:
\[
\begin{align}\frac{F'(x)}{F(x)}&=\frac{a_1+2a_2x+3a_3x^2+\cdots}{a_1x+a_2x^2+\cdots}\&=\frac{a_1+2a_2x+3a_3x^2+\cdots}{a_1x}\cdot \frac{1}{1+\left(\frac{a_2}{a_1}x+\frac{a_3}{a_1}x^2+\cdots\right)}\&=\left(x^{-1}+\frac{2a_2}{a_1}+\cdots\right)\left(1-x\left(\frac{a_2}{a_1}+\frac{a_3}{a_1}x+\cdots\right)\right)
\end{align}
\]
总之经过一系列代数变形得到 \([x^{-1}]\frac{F'(x)}{F(x)} = 1\)。
于是就有:
\[n*a_n =[x^{-1}](\frac{1}{F^n(x)})\]
\[a_n = \frac{1}{n}*[x^{-1}](\frac{1}{F^n(x)}) = \frac{1}{n}*[x^{n-1}](\frac{1}{(\frac{F(x)}{x})^n})\]
跑一个多项式逆元与多项式 exp + ln 求 k 次幂即可。
bzoj3684:大朋友和多叉树
题目点这里查看
根据题目给出的定义,设权值为 i 的结点为根的树方案数为 {ai},设 f(x) 为 {ai} 的生成函数。
则有:
\[f(x) = \sum_{k\in D}f^k(x) + x\]
其中最后一个 x 表示叶子结点的方案数。
稍微变形就可以得到复合逆的基本形式 \(g(f(x)) = x\),直接复合逆即可。
代码:
#include<cstdio>
#include<algorithm>
using namespace std;
const int MOD = 950009857;
const int MAXN = 400000;
const int G = 5;
int pow_mod(int b, int p) {
int ret = 1;
while( p ) {
if( p & 1 ) ret = 1LL*ret*b%MOD;
b = 1LL*b*b%MOD;
p >>= 1;
}
return ret;
}
struct poly{
int pw[20 + 5], ipw[20 + 5];
poly() {
for(int i=0;i<=21;i++)
pw[i] = pow_mod(G, (MOD-1)/(1<<i)), ipw[i] = pow_mod(pw[i], MOD-2);
}
void debug(int *A, int n) {
for(int i=0;i<n;i++)
printf("%d ", A[i]);
puts(""), puts("");
}
void clear(int *A, int l, int r) {
for(int i=l;i<r;i++)
A[i] = 0;
}
void copy(int *A, int *B, int n) {
for(int i=0;i<n;i++)
A[i] = B[i];
}
int length(int n) {
int len; for(len = 1; len < n; len <<= 1);
return len;
}
void ntt(int *A, int n, int type) {
for(int i=0,j=0;i<n;i++) {
if( i < j ) swap(A[i], A[j]);
for(int k=(n>>1);(j^=k)<k;k>>=1);
}
for(int i=1;(1<<i)<=n;i++) {
int s = (1<<i), t = (s>>1);
int u = (type == 1) ? pw[i] : ipw[i] ;
for(int j=0;j<n;j+=s) {
for(int k=0,p=1;k<t;k++,p=1LL*p*u%MOD) {
int x = A[j+k], y = 1LL*p*A[j+k+t]%MOD;
A[j+k] = (x + y)%MOD, A[j+k+t] = (x + MOD - y)%MOD;
}
}
}
if( type == -1 ) {
int iv = pow_mod(n, MOD-2);
for(int i=0;i<n;i++)
A[i] = 1LL*A[i]*iv%MOD;
}
}
int tmp1[MAXN + 5], tmp2[MAXN + 5];
void poly_mul(int *A, int *B, int *C, int n, int m, int k) {
int len = length(n + m - 1);
clear(tmp1, 0, len), clear(tmp2, 0, len);
for(int i=0;i<n;i++) tmp1[i] = A[i];
for(int i=0;i<m;i++) tmp2[i] = B[i];
ntt(tmp1, len, 1), ntt(tmp2, len, 1);
for(int i=0;i<len;i++)
C[i] = 1LL*tmp1[i]*tmp2[i]%MOD;
ntt(C, len, -1);
}
int tmp3[MAXN + 5];
void poly_inv(int *A, int *B, int n) {
if( n == 1 ) {
B[0] = pow_mod(A[0], MOD-2);
return ;
}
poly_inv(A, B, (n+1)>>1);
poly_mul(A, B, tmp3, n, (n+1)>>1, n);
for(int i=0;i<n;i++)
tmp3[i] = (MOD - tmp3[i])%MOD;
tmp3[0] = (tmp3[0] + 2)%MOD;
poly_mul(tmp3, B, B, n, (n+1)>>1, n);
}
void poly_int(int *A, int *B, int n) {
for(int i=n-1;i>0;i--)
B[i] = 1LL*pow_mod(i, MOD-2)*A[i-1]%MOD;
B[0] = 0;
}
void poly_der(int *A, int *B, int n) {
for(int i=0;i<n-1;i++)
B[i] = 1LL*(i+1)*A[i+1]%MOD;
B[n-1] = 0;
}
int tmp4[MAXN + 5];
void poly_ln(int *A, int *B, int n) {
poly_inv(A, tmp4, n);
poly_der(A, B, n);
poly_mul(B, tmp4, B, n, n, n);
poly_int(B, B, n);
}
int tmp5[MAXN + 5];
void poly_exp(int *A, int *B, int n) {
if( n == 1 ) {
B[0] = 1;
return ;
}
poly_exp(A, B, (n+1)>>1), clear(B, (n+1)>>1, n);
poly_ln(B, tmp5, n);
for(int i=0;i<n;i++)
tmp5[i] = (A[i] + MOD - tmp5[i])%MOD;
tmp5[0] = (tmp5[0] + 1)%MOD;
poly_mul(tmp5, B, B, n, (n+1)>>1, n);
}
int tmp6[MAXN + 5];
void poly_pow(int *A, int *B, int n, int k) {
poly_ln(A, tmp6, n);
for(int i=0;i<n;i++)
tmp6[i] = 1LL*k*tmp6[i]%MOD;
poly_exp(tmp6, B, n);
}
}op;
int s, m, f[MAXN + 5], g[MAXN + 5];
int main() {
scanf("%d%d", &s, &m);
for(int i=1;i<=m;i++) {
int x; scanf("%d", &x);
f[x] = (MOD - 1);
}
f[1] = 1;
for(int i=0;i<s;i++)
f[i] = f[i+1];
op.poly_pow(f, f, s, s);
op.poly_inv(f, g, s);
printf("%lld\n", 1LL*g[s-1]*pow_mod(s, MOD-2)%MOD);
}
形式化地给出问题:
若已知多项式函数 F(x),H(x),且复合函数 G(F(x)) 满足 G(F(x)) = H(x),求 G(x) 的第 n 项。
基本上和拉格朗日反演的推导过程是一样的,最后得到:
\[n*a_n =[x^{-1}](\frac{H'(x)}{F^n(x)})\]
\[a_n = \frac{1}{n}*[x^{n-1}](\frac{H'(x)}{(\frac{F(x)}{x})^n})\]
不过这里有一个常用的变换形式:
令 \(x' = F(x)\),则 \(x = F^{-1}(x')\)。这意味着我们可以将复合函数内部的某个函数“置换”出来。
咕。
标签:algorithm math mp3 += deb algo list printf 转换
原文地址:https://www.cnblogs.com/Tiw-Air-OAO/p/11233318.html