1. 程式人生 > 其它 >P5488-差分與字首和【NTT,生成函式】

P5488-差分與字首和【NTT,生成函式】

技術標籤:多項式luoguNTT生成函式

正題

題目連結:https://www.luogu.com.cn/problem/P5488


題目大意

求一個長度為 n n n的序列的 k k k階差分/字首和。


解題思路

先考慮字首和怎麼做
搞出來生成函式就是 ( ∑ i = 0 n a i x i ) ∗ ( ∑ i = 0 ∞ x i ) k (\sum_{i=0}^na_ix^i)*(\sum_{i=0}^{\infty}x^i)^k (i=0naixi)(i=0xi)k

然後根據常識我們知道 ( ∑ i = 0 ∞ x i ) k = ∑ i = 0 ∞ ( i + k − 1 i ) x i (\sum_{i=0}^{\infty}x^i)^k=\sum_{i=0}^{\infty}\binom{i+k-1}{i}x^i

(i=0xi)k=i=0(ii+k1)xi,當然也可以理解為 x i x^i xi的係數就是每次會往後跳任意格(可以是 0 0 0),然後 k k k次跳 i i i步的方案數。

之後式子 ( ∑ i = 0 n a i x i ) ∗ ( ∑ i = 0 ∞ ( i + k − 1 i ) x i ) (\sum_{i=0}^na_ix^i)*(\sum_{i=0}^{\infty}\binom{i+k-1}{i}x^i) (i=0naixi)(i=0(ii+k1)xi)

直接卷就好了

然後是差分,就是 ∑ i = 0 n a i x i ( ∑ i = 0 ∞ x i ) k \frac{\sum_{i=0}^na_ix^i}{(\sum_{i=0}^{\infty}x^i)^k}

(i=0xi)ki=0naixi
直接上多項式除法就好 。考慮怎麼轉換這個式子,我們需要用生成函式的方法了,上等比數列公式就有 ( ∑ i = 0 ∞ x i ) k = 1 ( 1 − x ) k (\sum_{i=0}^{\infty}x^i)^k=\frac{1}{(1-x)^k} (i=0xi)k=(1x)k1

那麼它的倒數就是 ( 1 − x ) k (1-x)^k (1x)k,上二項式定理就有 ∑ i = 0 k ( − 1 ) i ( k i ) x i \sum_{i=0}^k(-1)^i\binom{k}{i}x^i i=0k(1)i(ik)xi

然後也是兩個式子捲起來就好了。

然後 k k k可以直接模 p p p,考慮 f ( k ) = ( ∑ i = 0 ∞ x i ) k f(k)=(\sum_{i=0}^\infty x^i)^k f(k)=(i=0xi)k那麼有 f ( k n ) = f ( k ) n f(k^n)=f(k)^n f(kn)=f(k)n所以直接讓 k k k p p p即可

時間複雜度 O ( n log ⁡ n ) O(n\log n) O(nlogn)


c o d e code code

#include<cstdio>
#include<algorithm>
#include<cstring>
#define ll long long
using namespace std;
const ll N=(1e5)*6,P=1004535809,g=3;
struct poly{
    ll a[N],n;
}G,F;
ll n,k,t,inv[N],r[N];
char s[N];
ll power(ll x,ll b){
    ll ans=1;
    while(b){
        if(b&1)ans=ans*x%P;
        x=x*x%P;b>>=1;
    }
    return ans;
}
void NTT(ll *f,ll op,ll n){
    for(ll i=0;i<n;i++)
        if(r[i]<i)swap(f[i],f[r[i]]);
    for(ll p=2;p<=n;p<<=1){
        ll len=(p>>1);
        ll tmp=power(g,(P-1)/p);
        if(op)tmp=power(tmp,P-2);
        for(ll k=0;k<n;k+=p){
            ll buf=1;
            for(ll i=k;i<k+len;i++){
                ll tt=buf*f[i+len]%P;
                f[i+len]=(f[i]-tt+P)%P;
                f[i]=(f[i]+tt)%P;
                buf=buf*tmp%P;
            }
        }
    }
    if(op){
        int invn=power(n,P-2);
        for(int i=0;i<n;i++)
            f[i]=f[i]*invn%P;
    }
    return;
}
void mul(poly &a,poly &b){
    ll n=1;
    while(n<=a.n+b.n)n<<=1;
    for(ll i=0;i<n;i++)
        r[i]=(r[i>>1]>>1)^((i&1)?(n>>1):0);
    NTT(a.a,0,n);NTT(b.a,0,n);
    for(ll i=0;i<n;i++)
        a.a[i]=a.a[i]*b.a[i]%P;
    NTT(a.a,1,n);
    return;
}
int main()
{
    scanf("%lld",&n);
    scanf("%s",s);ll l=strlen(s);
    scanf("%lld",&t);
    for(ll i=0;i<l;i++)
        k=(k*10+s[i]-'0')%P;
    inv[1]=1;
    for(ll i=2;i<=n;i++)
        inv[i]=(P-(P/i)*inv[P%i]%P)%P;
    for(ll i=0;i<n;i++)
        scanf("%lld",&F.a[i]);
    G.a[0]=1;
    for(ll i=1;i<=n;i++){
        if(!t)G.a[i]=G.a[i-1]*(i+k-1)%P*inv[i]%P;
        else{
            if(i>k)break;
            G.a[i]=(-1)*G.a[i-1]*(k-i+1)%P*inv[i]%P;
            G.a[i]=(G.a[i]+P)%P;
        }
    }
    G.n=F.n=n;
    mul(G,F);
    for(ll i=0;i<n;i++)
        printf("%lld ",G.a[i]);
    return 0;
}