跳转至

习题

习题一

Problem - B - Codeforces

\[ A[i]=(\sum_{j=i}^{min(i+L-1,n)}A[j])\mod{P} \]

每次对 \(i\)\(1\)\(n\) 做上述式子,问做 \(m\) 次操作后 \(A[i]\) 的值?

\(A\) 翻转, 对于每次操作, 实际上就是与一个前 \(L\) 项都为 \(1\) 的多项式 \(f\) 进行卷积

我们可以使用多项式快速幂, 求的 \(f^m\), 最后计算 \(a * f^m\), 翻转输出即可

Code
C++
{ 多项式基础模版 }
void solve() {
    int n, L, m;
    cin >> n >> L >> m;
    Poly a(n);
    cin >> a;
    reverse(begin(a), end(a));
    Poly b(n);
    for (int i = 0 ; i < L ; i++) {
        b[i] = 1;
    }
    b = Pow(b, m);
    a *= b;
    a.resize(n);
    reverse(begin(a), end(a));
    cout << a << "\n";
}
int main() {
    int T;
    cin >> T;
    for (int i = 1 ; i <= T ; i++) {
        cout << "Case " << i << ": ";
        solve();
    }
    return 0;
}

P6115 【模板】整式递推 - 洛谷

求解递推系数是关于 \(n\) 的多项式的值。

\(M\) 是阶数。

\(D\) 是系数。

代码
C++
int M; // 阶数
int D; // 次数
struct Matrix {
    int a[10][10];
    Matrix() {
        memset(a, 0, sizeof a);
    }
    friend Matrix operator*(const Matrix &x, const Matrix &y) {
        Matrix res;
        for (int i = 0 ; i < M ; i++) {
            for (int j = 0 ; j < M ; j++) {
                for (int k = 0 ; k < M ; k++) {
                    res.a[i][j] += 1LL * x.a[i][k] * y.a[k][j] % P;
                    if (res.a[i][j] >= P) {
                        res.a[i][j] -= P;
                    }
                }
            }
        }
        return res;
    }
    friend bool operator==(const Matrix &x, const Matrix &y) {
        for (int i = 0 ; i < M ; i++) {
            for (int j = 0 ; j < M ; j++) {
                if (x.a[i][j] != y.a[i][j]) {
                    return false;
                }
            }
        }
        return true;
    }
};
Poly p[10];
Matrix get(i64 x) {
    x %= P;
    Matrix res;
    i64 p0 = eval(p[0], (x + M) % P);
    for (int i = 0 ; i < M - 1 ; i++) {
        res.a[i + 1][i] = p0;
    }
    for (int i = 0 ; i < M ; i++) {
        res.a[i][M - 1] = P - eval(p[M - i], (x + M) % P);
    }
    return res;
}
Matrix ff(int d, int x) {
    Matrix res = get(x);
    for (int i = 1 ; i < d;  i++) {
        res = res * get(x + i);
    }
    return res;
}
i64 hh(int d, int x) {
    i64 res = eval(p[0], x);
    for (int i = 1 ; i < d ; i++) {
        res = res * eval(p[0], x + i) % P;
    }
    return res;
}
const int N = 262145;
int fac[N], ifac[N];
void lagrange(const Matrix *F1, const int *F2, int n, i64 m, Matrix *R1, int *R2, bool flag) {
    static int pre[N], suf[N], f1[N], f2[N], g[N], inv[N], ifcm[N], facm[N];
    int k = n << 1 | 1, lim = norm((n << 1) + 1);
    facm[0] = 1;
    for (int i = 0 ; i <= n ; i++) {
        facm[0] = 1LL * facm[0] * (m - n + i) % P;
        ifcm[i] = 1LL * ifac[i] * ifac[n - i] % P;
    }
    pre[0] = suf[k + 1] = 1;
    for (int i = 1 ; i <= k ; i++) {
        pre[i] = 1LL * pre[i - 1] * (m - n + i - 1) % P;
    }
    for (int i = k ; i >= 1 ; i--) {
        suf[i] = 1LL * suf[i + 1] * (m - n + i - 1) % P;
    }
    i64 mul = qpow(pre[k], P - 2);
    for (int i = 1 ; i <= k ; i++) {
        inv[i] = mul * pre[i - 1] % P * suf[i + 1] % P;
    }
    for (int i = 1 ; i <= n ; i++) {
        facm[i] = 1LL * facm[i - 1] * (m + i) % P * inv[i] % P;
    }
    for (int i = 0 ; i < k ; i++) {
        g[i] = inv[i + 1];
    }
    for (int i = 0 ; i <= n ; i++) {
        f2[i] = 1LL * ifcm[i] * ((n - i) & 1 ? (P - F2[i]) : F2[i]) % P;
    }
    fill(g + k, g + lim + 1, 0);
    NTT::DIF(g, lim);
    fill(f2 + n + 1, f2 + lim + 1, 0);
    NTT::DIF(f2, lim);
    for (int i = 0 ; i < lim ; i++) {
        f2[i] = 1LL * f2[i] * g[i] % P;
    }
    NTT::IDIT(f2, lim);
    for (int i = 0 ; i <= n ; i++) {
        R2[i] = 1LL * f2[i + n] * facm[i] % P;;
    }
    if (!flag) {
        return;
    }
    for (int i = 0 ; i < M ; i++) {
        for (int j = 0 ; j < M ; j++) {
            for (int t = 0 ; t <= n ; t++) {
                f1[t] = 1LL * ifcm[t] * ((n - t) & 1 ? (P - F1[t].a[i][j]) : F1[t].a[i][j]) % P;
            }
            fill(f1 + n + 1, f1 + lim + 1, 0);
            NTT::DIF(f1, lim);
            for (int t = 0 ; t < lim ; t++) {
                f1[t] = 1LL * f1[t] * g[t] % P;
            }
            NTT::IDIT(f1, lim);
            for (int t = 0 ; t <= n ; t++) {
                R1[t].a[i][j] = 1LL * f1[t + n] * facm[t] % P;
            }
        }
    }
}
int st[31], top = 0;
int h[N], hd[N];
Matrix f[N], fd[N];
pair<Matrix, i64> magic(int s, int t) {
    int x = s;
    while(x) {
        st[++top] = x;
        x /= 2;
    }
    for (int i = 0 ; i <= D ; i++) {
        x = i * s;
        h[i] = eval(p[0], x);
        f[i] = get(x);
    }
    top--;
    int d = 1;
    int invs = qpow(s, P - 2);
    while (top--) {
        int Dd = D * d;
        lagrange(f, h, Dd, Dd + 1, f + Dd + 1, h + Dd + 1, true);
        h[Dd << 1 | 1] = 0;
        f[Dd << 1 | 1] = Matrix();
        lagrange(f, h, Dd << 1, 1LL * d * invs % P, fd, hd, true);
        for (int i = 0 ; i <= (Dd << 1) ; i++) {
            f[i] = f[i] * fd[i];
            h[i] = 1LL * h[i] * hd[i] % P;
        }
        d <<= 1;
        if (!(st[top + 1] & 1)) {
            continue;
        }
        Dd = D * (d + 1);
        for (int i = D * d + 1 ; i <= Dd ; i++) {
            x = i * s;
            f[i] = ff(d, x);
            h[i] = hh(d, x);
        }
        for (int i = 0 ; i <= Dd ; i++) {
            x = i * s;
            f[i] = f[i] * get(x + d);
            h[i] = h[i] * eval(p[0], x + d) % P;
        }
        d++;
    }
    lagrange(f, h, 1LL * D * s, 1LL * M * invs % P, f, h, false);
    Matrix r1;
    for (int i = 0 ; i < M ; i++) {
        r1.a[i][i] = 1;
    }
    i64 r2 = 1;
    for (int i = 0 ; i <= t ; i++) {
        r1 = r1 * f[i];
        r2 = r2 * h[i] % P;
    }
    return {r1, r2};
}
i64 recursive(const int *a, i64 n) {
    i64 tn = n - M + 1;
    i64 s = ceill(sqrtl(tn * 1.0 / D)) + 1;
    i64 res = 0;
    auto [mul, div] = magic(s, (tn - s) / s);
    for (i64 i = (tn / s) * s ; i < tn ; i++) {
        mul = mul * get(i % P);
        div = div * eval(p[0], i + M) % P;
    }
    for (int i = 0 ; i < M ; i++) {
        res = (res + 1LL * a[i] * mul.a[i][M - 1] % P) % P;
    }
    return res * qpow(div, P - 2) % P;
}
int a[10];
int main() {
    i64 n;
    int m, d;
    cin >> n >> m >> d;
    for (int i = 0 ; i < m ; i++) {
        cin >> a[i];
    }
    for (int i = 0 ; i <= m ; i++) {
        p[i].resize(d + 1);
        for (int j = 0 ; j <= d ; j++) {
            cin >> p[i][j];
        }
    }
    {
        // 整式递推初始化
        // M 是阶数,对应 a 数组, 表示前 M 项的值, 
        // D 是系数多项式的最高次,对应 p 数组, 表示系数多项式.
        M = m;
        D = d;
        fac[0] = 1;
        for (int i = 1 ; i < N ; i++) {
            fac[i] = 1LL * fac[i - 1] * i % P;
        }
        ifac[N - 1] = qpow(fac[N - 1], P - 2);
        for (int i = N - 1 ; i >= 0 ; i--) {
            ifac[i - 1] = 1LL * ifac[i] * i % P;
        }
    }
    i64 ans = recursive(a, n);
    cout << ans;
    return 0;
}