• [多项式][洛谷 P5488] 差分与前缀和


    题目大意

    给定一个长为 (n) 的序列 ({a_n}),求出其 (k) 阶差分或前缀和序列。
    结果的每一项都需要对 (1004535809) 取模。(1leq nleq 10^5,0leq a_ileq 10^9,1leq kleq 10^{2333},k otequiv 0(mathrm{mod} 1004535809))

    题解

    首先考虑前缀和。

    ({a_n}) 的一阶前缀和序列为 ({s_n}),设 ({a_n}) 的生成函数 (a(x)=sum_{i=1}^n a_ix^i)({s_n}) 的生成函数是 (s(x)=sum_{i=1}^n s_ix^i)

    (s_n=x^nsum_{i=1}^na_i=sum_{i=1}^na_ix^icdot x^{n-i+1})

    (g(x)=sum_{i=0}^{infty} x^i),则 (s(x)=a(x)*g(x))

    那么 ({a_n})(k) 阶前缀的生成函数(s_k(x)=a(x)*g^k(x))

    由麦克劳林公式, (sum_{i=0}^{infty}x^i=frac{1}{1-x}),所以 (s_k(x)=a(x)*left(frac{1}{1-x} ight)^k)

    考虑怎么求 (left(frac{1}{1-x} ight)^k)(n) 次项系数。由广义二项式定理有

    [left(frac{1}{1-x} ight)^k=(1-x)^{-k}=sum_{i=0}^{infty}inom{-k}{i}(-x)^i\ =sum_{i=0}^{infty}frac{-k^{underline{i}}}{i!}(-x)^i=sum_{i=0}^{infty}frac{(k+i-1)^{underline{i}}}{i!}x^i\ =sum_{i=0}^{infty}inom{k+i-1}{i}x^i ]

    所以 (left(frac{1}{1-x} ight)^k)(n) 次项系数为 (inom{k+n-1}{n})

    (f(n)=inom{k+n-1}{n}),则 (f(n)=frac{k+n-1}{n}f(n-1))
    我们先求出第一项,然后递推求出后面的系数。

    然后使用NTT对 (a(x))(left(frac{1}{1-x} ight)^k) 作卷积,即可求出 (k) 阶前缀和。

    然后考虑差分。

    (h(x))({a_n}) 的一阶差分的生成函数,则

    [h(x)=sum_{i=1}^n(a_i-a_{i-1})x^iequivsum_{i=1}^na_ix^i-xsum_{i=1}^na_ix^iequiv(1-x)*a(x) pmod{x^{n+1}} ]

    那么 (k) 阶差分有

    [h_k(x)=a(x)*(1-x)^k ]

    由二项式定理

    [(1-x)^k=sum_{i=1}^kinom{k}{i}(-x)^k=sum_{i=1}^k(-1)^kinom{k}{i} x^i ]

    于是 ((1-x)^k)(n) 次项系数是 ((-1)^kinom{k}{i})

    (inom{k}{i}=frac{k-n+1}{n}inom{k}{n-1}),所以我们也可以递推求得。

    Code

    #include <bits/stdc++.h>
    using namespace std;
    
    #define RG register int
    #define LL long long
    
    template<typename elemType>
    inline void Read(elemType& T) {
        elemType X = 0, w = 0; char ch = 0;
        while (!isdigit(ch)) { w |= ch == '-';ch = getchar(); }
        while (isdigit(ch)) X = (X << 3) + (X << 1) + (ch ^ 48), ch = getchar();
        T = (w ? -X : X);
    }
    
    inline LL ExPow(LL b, LL n, LL MOD) {
        if (MOD == 1) return 0;
        LL x = 1;
        LL Power = b % MOD;
        while (n) {
            if (n & 1) x = x * Power % MOD;
            Power = Power * Power % MOD;
            n >>= 1;
        }
        return x;
    }
    
    const int maxn = 270000;
    const LL P = 1004535809LL, G = 3, Gi = 334845270;
    
    namespace Poly {
        int r[maxn];
        LL buf[maxn], buf1[maxn], buf2[maxn], buf3[maxn], buf4[maxn];
        int L, limit;
    
        LL pinv(LL x) { return ExPow(x, P - 2, P); }
        //快速数论变换 type=1:正变换 type=-1:逆变换
        void NTT(LL* A, int type) {
            for (int i = 0; i < limit; i++)
                if (i < r[i]) swap(A[i], A[r[i]]);
            for (int mid = 1; mid < limit; mid <<= 1) {
                LL Wn = ExPow(type == 1 ? G : Gi, (P - 1) / (mid << 1), P);
                for (int j = 0; j < limit; j += (mid << 1)) {
                    LL w = 1;
                    for (int k = 0; k < mid; k++, w = (w * Wn) % P) {
                        int x = A[j + k], y = w * A[j + k + mid] % P;
                        A[j + k] = (x + y) % P;
                        A[j + k + mid] = (x - y + P) % P;
                    }
                }
            }
            if (type == 1) return;
            LL inv_limit = pinv(limit);
            for (int i = 0; i < limit; ++i)
                A[i] = A[i] * inv_limit % P;
        }
    
        //多项式卷积 a(x): N-1次多项式 b(x): M-1次多项式
        void Conv(LL* A, int N, LL* B, LL M, LL* C) {
            L = 0; limit = 1;
            while (limit <= N + M) limit <<= 1, L++;
            for (int i = 0; i < limit; i++) r[i] = (r[i >> 1] >> 1) | ((i & 1) << (L - 1));
            NTT(A, 1); NTT(B, 1);
            for (int i = 0; i < limit; i++) C[i] = A[i] * B[i] % P;
            NTT(C, -1);
        }
    }
    
    char K[2500];
    LL a[maxn], b[maxn], c[maxn], s[maxn], h[maxn], Inv[maxn];
    int n, m, opt;
    LL k = 0;
    
    void Init() {
        Inv[1] = 1;
        for (int i = 2;i <= 200000;++i)
            Inv[i] = ((-(P / i) * Inv[P % i]) % P + P) % P;
        return;
    }
    
    void PrefixSum() { // k阶前缀和
        s[0] = 1;
        for (int i = 1;i <= n;++i)
            s[i] = ((k + i - 1) % P + P) % P * Inv[i] % P * s[i - 1] % P;
        Poly::Conv(a, n + 1, s, n + 1, s);
        for (int i = 1;i <= n;++i)
            printf("%lld ", s[i]);
        printf("
    ");
    }
    
    void Diff() { // k阶差分
        h[0] = 1;
        for (int i = 1;i <= n;++i)
            h[i] = (-1LL * (k - i + 1) % P + P) % P * Inv[i] % P * h[i - 1] % P;
        Poly::Conv(a, n + 1, h, n + 1, h);
        for (int i = 1;i <= n;++i)
            printf("%lld ", h[i]);
        printf("
    ");
    }   
    
    int main() {
        Init();
        scanf("%d%s%d", &n, K + 1, &opt);
        m = strlen(K + 1);
        for (int i = 1;i <= m;++i)
            k = (k * 10 + K[i] - '0') % P;
        for (int i = 1;i <= n;++i)
            Read(a[i]);
        if (opt == 0) PrefixSum();
        else Diff();
            
        return 0;
    }
    
  • 相关阅读:
    Win7中隐藏的上帝模式——GodMode
    C# 中有关 using 关键字
    数据结构 实验三  二叉树
    数据结构 实验二 栈
    指针知识(六):指针的指针
    指针知识(五):指针(pointer)和常量(const)
    指针知识(四):指针数学计算
    指针知识(三):指针与数组
    指针知识(二):指针初始化
    指针知识(一):指针声明
  • 原文地址:https://www.cnblogs.com/AEMShana/p/14563800.html
Copyright © 2020-2023  润新知