【luogu P6800】Chirp Z-Transform(多项式)(NTT)(bluestein 算法)

Chirp Z-Transform

题目链接:luogu P6800

题目大意

给你一个多项式和 c,m,要你求把 c^0,c^1,...c,^m-1 分别带入多项式得到的值。

思路

考虑把答案也看做是多项式:
ansi=F(ci)=j=0n1ajcij

然后又一个东西是:ij=(i+j2)(i2)(j2)
简单证明:
(i+j2)(i2)(j2)=(i+j)(i+j1)i(i1)j(j1)2
=i2+iji+ij+j2ji2+ij2+j2=2ij2=ij

然后带入:
ansi=F(ci)=j=0n1ajcij
=j=0n1ajc(i+j2)(i2)(j2)
=c(i2)j=0n1ajc(i+j2)c(j2)

然后发现右边这个部分可以直接卷起来,可以用 NTT 搞。

然后接着是求 c 的要用的次方项,每次都 log 太慢了,我们可以用预处理光速乘或者两边前缀和搞得出它。
两边阶乘的原理是 x(x1)2=(x1)((x1)+1)2=1+2+...+n1,然后就可以两边前缀和搞了。

你可以 ajc(j2) 弄一个多项式(系数变成 nj),然后 c(i+j2) 弄一个多项式(系数就是这个),然后加了之后就是我们要的第 i 项了。(外面那个 c(i2) 最后输出的时候乘上即可)

代码

#include<cstdio> #include<iostream> #include<algorithm> #define ll long long #define mo 998244353 using namespace std; int n, m, nm, limit, l_size; int an[3000001]; ll a[3000001], c, invc, G, Gv; ll jc[3000001], inv[3000001]; ll g[3000001]; ll ksm(ll x, ll y) { ll re = 1; while (y) { if (y & 1) re = re * x % mo; x = x * x % mo; y >>= 1; } return re; } void NTT(ll *now, int op) {//NTT for (int i = 0; i < limit; i++) if (i < an[i]) swap(now[i], now[an[i]]); for (int mid = 1; mid < limit; mid <<= 1) { ll Wn = ksm(op == 1 ? G : Gv, (mo - 1) / (mid << 1)); for (int R = (mid << 1), j = 0; j < limit; j += R) { ll w = 1; for (int k = 0; k < mid; k++, w = w * Wn % mo) { ll x = now[j + k], y = w * now[j + mid + k] % mo; now[j + k] = (x + y) % mo; now[j + mid + k] = (x - y + mo) % mo; } } } } int main() { scanf("%d %lld %d", &n, &c, &m); nm = max(n, m); invc = ksm(c, mo - 2);//两次前缀和求出 n*(n-1)/2 次方的阶乘以及它的逆元 jc[0] = 1; for (int i = 1; i <= n + m; i++) jc[i] = jc[i - 1] * c % mo; for (int i = 1; i <= n + m; i++) jc[i] = jc[i - 1] * jc[i] % mo; inv[0] = 1; for (int i = 1; i <= nm; i++) inv[i] = inv[i - 1] * invc % mo; for (int i = 1; i <= nm; i++) inv[i] = inv[i - 1] * inv[i] % mo; for (int i = 0; i < n; i++) scanf("%lld", &a[n - i]), a[n - i] = a[n - i] * (i ? inv[i - 1] : 1) % mo; for (int i = 0; i <= n + m; i++) g[i] = i ? jc[i - 1] : 1; limit = 1; while (limit <= n + m) { limit <<= 1; l_size++; } for (int i = 0; i < limit; i++) an[i] = (an[i >> 1] >> 1) | ((i & 1) << (l_size - 1)); G = 3; Gv = ksm(G, mo - 2); NTT(a, 1); NTT(g, 1); for (int i = 0; i < limit; i++) a[i] = a[i] * g[i] % mo; NTT(a, -1); ll liv = ksm(limit, mo - 2); for (int i = 0; i < limit; i++) a[i] = a[i] * liv % mo; for (int i = 0; i < m; i++) { printf("%lld ", a[i + n] * (i ? inv[i - 1] : 1) % mo); } return 0; }

__EOF__

本文作者あおいSakura
本文链接https://www.cnblogs.com/Sakura-TJH/p/luogu_P6800.html
关于博主:评论和私信会在第一时间回复。或者直接私信我。
版权声明:本博客所有文章除特别声明外,均采用 BY-NC-SA 许可协议。转载请注明出处!
声援博主:如果您觉得文章对您有帮助,可以点击文章右下角推荐一下。您的鼓励是博主的最大动力!
posted @   あおいSakura  阅读(49)  评论(0编辑  收藏  举报
编辑推荐:
· go语言实现终端里的倒计时
· 如何编写易于单元测试的代码
· 10年+ .NET Coder 心语,封装的思维:从隐藏、稳定开始理解其本质意义
· .NET Core 中如何实现缓存的预热?
· 从 HTTP 原因短语缺失研究 HTTP/2 和 HTTP/3 的设计差异
阅读排行:
· 分享一个免费、快速、无限量使用的满血 DeepSeek R1 模型,支持深度思考和联网搜索!
· 基于 Docker 搭建 FRP 内网穿透开源项目(很简单哒)
· ollama系列01:轻松3步本地部署deepseek,普通电脑可用
· 25岁的心里话
· 按钮权限的设计及实现
点击右上角即可分享
微信分享提示