P5205 【模板】多项式开根

\(\color{#0066ff}{ 题目描述 }\)

给定一个\(n-1\)次多项式\(A(x)\)，求一个在\(\bmod\ x^n\)意义下的多项式\(B(x)\)，使得\(B^2(x) \equiv A(x) \ (\bmod\ x^n)\)

多项式的系数在\(\bmod\ 998244353\)的意义下进行运算。

\(\color{#0066ff}{输入格式}\)

第一行一个正整数\(n\)。

接下来\(n\)个整数，依次表示多项式的系数\(a_0, a_1, \dots, a_{n-1}\)

保证\(a_0 = 1\).

\(\color{#0066ff}{输出格式}\)

输出\(n\)个整数，表示答案多项式的系数\(b_0, b_1, \dots, b_{n-1}\)

\(\color{#0066ff}{输入样例}\)

3
1 2 1

    
7
1 8596489 489489 4894 1564 489 35789489  

\(\color{#0066ff}{输出样例}\)

1 1 0

    
1 503420421 924499237 13354513 217017417 707895465 411020414

\(\color{#0066ff}{数据范围与提示}\)

对于\(100\%\)的数据：\(n \leq 10^5 \qquad a_i \in [0,998244352] \cap \mathbb{Z}\)

\(\color{#0066ff}{题解}\)

牛顿迭代走一波（^_~）

\(F(x)\equiv \sqrt {A(x)}\)

\(G(F(x))=F^2(x)-A(x)\)

\(G'(F(x))=2F(x)\)

\(F(x)\equiv F_0(x)-\frac{F_0^2(x)-A(x)}{2F_0(x)}=\frac{F_0^2(x)+A(x)}{2F_0(x)}\)

直接多项式求逆递归即可

// luogu-judger-enable-o2
#include<bits/stdc++.h>
#define LL long long
LL in() {
    char ch; LL x = 0, f = 1;
    while(!isdigit(ch = getchar()))(ch == '-') && (f = -f);
    for(x = ch ^ 48; isdigit(ch = getchar()); x = (x << 1) + (x << 3) + (ch ^ 48));
    return x * f;
}
const int maxn = 4e6 + 10;
const int mod = 998244353;
using std::vector;
int r[maxn], len;
LL ksm(LL x, LL y) {
    LL re = 1LL;
    while(y) {
        if(y & 1) re = re * x % mod;
        x = x * x % mod;
        y >>= 1;
    }
    return re;
}
void FFT(vector<int> &A, int flag) {
    A.resize(len);
    for(int i = 0; i < len; i++) if(i < r[i]) std::swap(A[i], A[r[i]]);
    for(int l = 1; l < len; l <<= 1) {
        int w0 = ksm(3, (mod - 1) / (l << 1));
        for(int i = 0; i < len; i += (l << 1)) {
            int w = 1, a0 = i, a1 = i + l;
            for(int k = 0; k < l; k++, a0++, a1++, w = 1LL * w0 * w % mod) {
                int tmp = 1LL * A[a1] * w % mod;
                A[a1] = ((A[a0] - tmp) % mod + mod) % mod;
                A[a0] = (A[a0] + tmp) % mod;
            }
        }
    }
    if(!(~flag)) {
        int inv = ksm(len, mod - 2);
        std::reverse(A.begin() + 1, A.end());
        for(int i = 0; i < len; i++) A[i] = 1LL * A[i] * inv % mod;
    }
}
vector<int> operator * (vector<int> A, vector<int> B) {
    int tot = A.size() + B.size() - 1;
    for(len = 1; len <= tot; len <<= 1);
    for(int i = 0; i < len; i++) r[i] = (r[i >> 1] >> 1) | ((i & 1) * (len >> 1));
    FFT(A, 1), FFT(B, 1);
    vector<int> ans;
    for(int i = 0; i < len; i++) ans.push_back(1LL * A[i] * B[i] % mod);
    FFT(ans, -1);
    ans.resize(tot);
    return ans;
}
vector<int> operator + (const vector<int> &A, const vector<int> &B) {
    vector<int> ans;
    for(int i = 0; i < (int)std::min(A.size(), B.size()); i++) ans.push_back(A[i] + B[i]);
    for(int i = A.size(); i < (int)B.size(); i++) ans.push_back(B[i]);
    for(int i = B.size(); i < (int)A.size(); i++) ans.push_back(A[i]);
    return ans;
}
vector<int> operator - (const vector<int> &A, const vector<int> &B) {
    vector<int> ans;
    for(int i = 0; i < (int)std::min(A.size(), B.size()); i++) ans.push_back(A[i] - B[i]);
    for(int i = A.size(); i < (int)B.size(); i++) ans.push_back(-B[i]);
    for(int i = B.size(); i < (int)A.size(); i++) ans.push_back(A[i]);
    return ans;
}
vector<int> inv(const vector<int> &A) {
    if(A.size() == 1) {
        vector<int> ans;
        ans.push_back(ksm(A[0], mod - 2));
        return ans;
    }
    vector<int> ans, B = A;
    int n = A.size(), _ = (n + 1) >> 1;
    B.resize(_);
    B = inv(B);
    ans.push_back(2);
    ans = B * (ans - A * B);
    ans.resize(n);
    return ans;
}
vector<int> sqt(const vector<int> &A) {
    if(A.size() == 1) return A;
    vector<int> ans, B = A;
    int n = A.size(), _ = (n + 1) >> 1;
    B.resize(_);
    B = sqt(B);
    B.resize(n, 0);
    ans.push_back(2);
    ans = (B * B + A) * inv(ans * B);
    ans.resize(n);
    return ans;
}
        
int main() {
    int n = in();
    vector<int> a;
    for(int i = 1; i <= n; i++) a.push_back(in());
    a = sqt(a);
    for(int i = 0; i < n; i++) printf("%d%c", a[i], i == n - 1? '\n' : ' ');
    return 0;
}