Loading [MathJax]/jax/output/CommonHTML/fonts/TeX/Typewriter-Regular.js

一些闲话 & 免责声明:感觉好久没有更过学习笔记了,现在才来学这个不知道来不来得及呢。因为笔者是一个数学菜鸡,所以并不会很严谨(。

问题描述

Portal.Portal.

LSB-first 算法

从线性递推式到分式第 nn

Q(x)=1f1xf2x2fkxkQ(x)=1f1xf2x2fkxkF(x)F(x) 为递推数列的 OGF。我们需要求解 F(x)n 次项系数。这里直接给出一个结论:

存在一个最高次项的次数小于 k 的多项式 P(x) 满足

F(x)=P(x)Q(x)

由于 P(x) 的度数 <k,而 F(x) 的前 k 项是已知的,所以可以仅用 F(x) 的前 k 项与 Q(x) 进行多项式乘法求出 P(x).

剩下的问题是,要怎么快速计算 P(x)Q(x) 的第 n 项。

Bostan-mori 算法

首先有

[xn]P(x)Q(x)=[xn]P(x)Q(x)Q(x)Q(x)

可以发现,Q(x)Q(x) 是一个偶函数(考虑将 x,x 分别代入),也就是说,它的奇次项均为零,于是可以构造多项式 v(x2)=Q(x)Q(x),同时对分子奇偶分类可得

[xn]ceven(x2)+xcodd(x2)v(x2)={[xn/2]ceven(x)v(x),if n is even[xn/2]codd(x)v(x),if n is odd

对于后面的一步推导,考虑 v(x) 的逆也是偶函数即可。于是可以进行递推求解。复杂度 O(klogklogn).

代码实现

三个 998244353 相加又忘记开 long long 了,我自裁 🔫。

/* 先复习动人的墓志,再背诵浪漫的情诗。 */

# include <cstdio>
# include <cctype>
# define print(x,y) write(x), putchar(y)

template <class T>
inline T read(const T sample) {
    T x=0; char s; bool f=0;
    while(!isdigit(s=getchar())) f|=(s=='-');
    for(; isdigit(s); s=getchar()) x=(x<<1)+(x<<3)+(s^48);
    return f? -x: x;
}
template <class T>
inline void write(T x) {
    static int writ[50], w_tp=0;
    if(x<0) putchar('-'), x=-x;
    do writ[++w_tp]=x-x/10*10, x/=10; while(x);
    while(putchar(writ[w_tp--]^48), w_tp);
}

# include <cstring>
# include <iostream>
using namespace std;

const int MAXN = 170000;
const int mod = 998244353, G = 3;

void mul(int& x,int y) { x=1ll*x*y%mod; }
int inc(int x,int y) { return x+y>=mod?x+y-mod:x+y; }
int dec(int x,int y) { return x-y<0?x-y+mod:x-y; }
int inv(int x,int y=mod-2,int r=1) {
    for(; y; y>>=1, x=1ll*x*x%mod)
        if(y&1) r=1ll*r*x%mod; return r;
} const int iG = inv(G);

int P[MAXN], rev[MAXN], g[MAXN];
int lim, bit, Inv, n, k, Q[MAXN];

void ntt(int* f,bool opt=0) {
    for(int i=0;i<lim;++i)
        if(i<rev[i]) swap(f[i],f[rev[i]]);
    for(int mid=1;mid<lim;mid<<=1) {
        const int wn = inv(opt?iG:G,(mod-1)/(mid<<1));
        for(int i=0; i<lim; i+=(mid<<1)) { int w=1;
            for(int j=0; j<mid; ++j, w=1ll*w*wn%mod) {
                const int tmp = 1ll*f[i|j|mid]*w%mod;
                f[i|j|mid] = dec(f[i|j], tmp),
                f[i|j] = inc(f[i|j], tmp);
            }
        }
    } if(!opt) return;
    for(int i=0;i<lim;++i) mul(f[i],Inv);
}

int main() {
    n=read(9), k=read(9); Q[0]=lim=1;
    for(int i=1;i<=k;++i) Q[i] = (0ll+(mod<<1)-read(9))%mod;
    for(int i=0;i<k;++i) P[i] = (0ll+(mod<<1)+read(9))%mod;
    while(lim<(k<<1|1)) lim<<=1, ++bit; 
    for(int i=0;i<lim;++i)
        rev[i] = (rev[i>>1]>>1)|((i&1)<<(bit-1)); 
    ntt(P), ntt(Q); Inv = inv(lim);
    for(int i=0;i<lim;++i) mul(P[i],Q[i]);
    ntt(P,1), ntt(Q,1);
    const int SIZE = sizeof(int)*(lim-k-1);
    memset(P+k,0,sizeof(int)*(lim-k));
    memset(Q+k+1,0,SIZE);
    for(; n; n>>=1) {
        for(int i=0;i<=k;++i)
            g[i] = (i&1)?(mod-Q[i])%mod:Q[i];
        memset(g+k+1,0,SIZE);
        ntt(g), ntt(P), ntt(Q);
        for(int i=0;i<lim;++i)
            mul(P[i],g[i]), mul(Q[i],g[i]);
        ntt(P,1), ntt(Q,1);
        for(int i=0;i<=k;++i) Q[i]=Q[i<<1];
        for(int i=0;i<=k;++i) P[i]=P[i<<1|(n&1)];
        memset(P+k+1,0,SIZE), memset(Q+k+1,0,SIZE);
    } print((mul(P[0],inv(Q[0])), P[0]), '\n');
    return 0;                               
}
posted on   Oxide  阅读(241)  评论(0编辑  收藏  举报
编辑推荐:
· DeepSeek 解答了困扰我五年的技术问题
· 为什么说在企业级应用开发中,后端往往是效率杀手?
· 用 C# 插值字符串处理器写一个 sscanf
· Java 中堆内存和栈内存上的数据分布和特点
· 开发中对象命名的一点思考
阅读排行:
· DeepSeek 解答了困扰我五年的技术问题。时代确实变了!
· PPT革命!DeepSeek+Kimi=N小时工作5分钟完成?
· What?废柴, 还在本地部署DeepSeek吗?Are you kidding?
· DeepSeek企业级部署实战指南:从服务器选型到Dify私有化落地
· 程序员转型AI:行业分析
历史上的今天:
2021-07-22 [杂题合集] One Last Blog



点击右上角即可分享
微信分享提示