OI

Universal OJ#62. 「UR #5」怎样跑得更快 – 题解

主要思路

优美的神仙操作,爱了爱了。

把原式抄下来:

\[ \begin{equation} \sum_{j = 1}^{n} \gcd(i, j)^c \cdot \text{lcm}(i, j)^d \cdot x_j \equiv b_i \pmod{p} \end{equation} \]

首先把烦人的 \(\text{lcm}\) 去掉:

\[ \begin{aligned} \sum_{j = 1}^n \gcd(i, j)^{c – d} i^d j^d \cdot x_j &\equiv b_i \pmod p \\ \sum_{j = 1}^n \gcd(i, j)^{c – d} j^d \cdot x_j &\equiv \frac{b_i}{i^d} \pmod p \end{aligned} \]

再进行一些标记。设 \(X_i = x_i i^d, B_i = \frac{b_i}{i^d}\)。得到:

\[ \sum_{j = 1}^n \gcd(i, j)^{c – d} X_j \equiv B_i \pmod p \]

考虑把这个烦人的 \(\gcd\) 给弄掉。考虑经典的莫比乌斯反演:

\[ \begin{aligned} f(x) &= x^{c – d} = \sum_{k|x} f'(k) \\ f'(x) &= \sum_{k|x} \mu(\frac{x}{k}) f(k) \\ \sum_{j = 1}^n X_j f'(\gcd(i, j)) &\equiv B_i \pmod p \end{aligned} \]

试试其他枚举方式:

\[ \begin{aligned} \sum_{j = 1}^n X_j \sum_{d|i, d|j} f'(d) &\equiv B_i \pmod p \\ \sum_{d|i} f'(d) \sum_{d|j} X_j &\equiv B_i \pmod p \end{aligned} \]

设 \(h(d) = \sum_{d|j} X_j\),得:

\[ \begin{aligned} \sum_{d|i} f'(d) h(d) &\equiv B_i \pmod p \\ f'(i)h(i) &\equiv B_i – \sum_{d|i, d \neq i} f'(d)h(d) \end{aligned} \]

调和级数一波即可。求得 \(f'(d)h(d)\) 之后先去掉 \(f'(d)\),再来分解:

\[ \begin{aligned} h(d) &= \sum_{d|j} X_j \\ X_d &= h(d) – \sum_{d|j, j \neq d} X_j \end{aligned} \]

基本差不多了。无解的情况只需要中间判一判即可。

代码

// UOJ62.cpp
#include <bits/stdc++.h>

using namespace std;

#define int long long

const int MAX_N = 1e5 + 200, mod = 998244353;

int n, c, d, q, bi[MAX_N], mu[MAX_N], primes[MAX_N], ptot, pc[MAX_N], pd[MAX_N], fd[MAX_N], h[MAX_N], xi[MAX_N];
bool vis[MAX_N];

int fpow(int bas, int tim)
{
    if (bas == 0 && tim == 0)
        return 0;
    int ret = 1;
    while (tim)
    {
        if (tim & 1)
            ret = 1LL * ret * bas % mod;
        bas = 1LL * bas * bas % mod;
        tim >>= 1;
    }
    return ret;
}

void sieve()
{
    mu[1] = 1;
    for (int i = 2; i < MAX_N; i++)
    {
        if (!vis[i])
            primes[++ptot] = i, mu[i] = mod - 1;
        for (int j = 1; j <= ptot && 1LL * i * primes[j] < MAX_N; j++)
        {
            vis[i * primes[j]] = true;
            if (i % primes[j] == 0)
                break;
            mu[i * primes[j]] = mod - mu[i];
        }
    }
    for (int i = 1; i < MAX_N; i++)
        pc[i] = fpow(i, c), pd[i] = fpow(i, d);
    for (int i = 1; i < MAX_N; i++)
    {
        int f = 1LL * pc[i] * fpow(pd[i], mod - 2) % mod;
        for (int j = 1; i * j < MAX_N; j++)
            fd[i * j] = (0LL + fd[i * j] + 1LL * mu[j] * f % mod) % mod;
    }
}

signed main()
{
    scanf("%lld%lld%lld%lld", &n, &c, &d, &q), sieve();
    while (q--)
    {
        for (int i = 1; i <= n; i++)
            scanf("%lld", &bi[i]), bi[i] = 1LL * bi[i] * fpow(pd[i], mod - 2) % mod;
        bool flag = true;
        for (int i = 1; i <= n; i++)
        {
            if (fd[i] == 0 && bi[i] != 0)
            {
                flag = false;
                break;
            }
            h[i] = 1LL * fpow(fd[i], mod - 2) * bi[i] % mod;
            for (int j = 2 * i; j <= n; j += i)
                bi[j] = (0LL + bi[j] + mod - bi[i]) % mod;
        }
        if (!flag)
        {
            puts("-1");
            continue;
        }
        for (int i = n; i; i--)
        {
            xi[i] = h[i];
            for (int j = 2 * i; j <= n; j += i)
                xi[i] = (0LL + xi[i] + mod - xi[j]) % mod;
        }
        for (int i = 1; i <= n; i++)
            printf("%lld ", 1LL * xi[i] * fpow(pd[i], mod - 2) % mod);
        puts("");
    }
    return 0;
}

 


kal0rona

http://kaloronahuang.com

江西师大附中全机房最弱