P4213 【模板】杜教筛、P3768 简单的数学题、P3803 【模板】多项式乘法（FFT）

P4213 【模板】杜教筛

题目描述

P4213 【模板】杜教筛 - 洛谷 | 计算机科学教育新生态 (luogu.com.cn)

运行代码

#include <iostream>
#include <map>

using namespace std;

const int N = 2000010;
long long vis[N], pm[N], mu[N], phi[N], cnt;
map<long long, long long> mp_mu, mp_phi;

void init() {
    mu[1] = phi[1] = 1;
    for (int i = 2; i < N; i++) {
        if (!vis[i]) {
            pm[++cnt] = i;
            mu[i] = -1;
            phi[i] = i - 1;
        }
        for (int j = 1; i * pm[j] < N; j++) {
            int p = pm[j];
            vis[i * p] = 1;
            if (i % p == 0) {
                phi[i * p] = phi[i] * p;
                break;
            }
            mu[i * p] = -mu[i];
            phi[i * p] = phi[i] * (p - 1);
        }
    }
    for (int i = 1; i < N; i++) {
        mu[i] += mu[i - 1];
        phi[i] += phi[i - 1];
    }
}

long long Sphi(long long n) {
    if (n < N) {
        return phi[n];
    }
    if (mp_phi.count(n)) {
        return mp_phi[n];
    }
    long long ans = n * (n + 1) / 2;
    for (long long l = 2, r; l <= n; l = r + 1) {
        r = n / (n / l);
        ans -= Sphi(n / l) * (r - l + 1);
    }
    return mp_phi[n] = ans;
}

long long Smu(long long n) {
    if (n < N) {
        return mu[n];
    }
    if (mp_mu.count(n)) {
        return mp_mu[n];
    }
    long long ans = 1;
    for (long long l = 2, r; l <= n; l = r + 1) {
        r = n / (n / l);
        ans -= Smu(n / l) * (r - l + 1);
    }
    return mp_mu[n] = ans;
}

int main() {
    init();
    long long T, n;
    cin >> T;
    while (T--) {
        cin >> n;
        cout << Sphi(n) << " " << Smu(n) << endl;
    }
    return 0;
}

代码思路

函数 init ：首先对 mu[1] 和 phi[1] 进行初始化，将它们都设为 1 。然后通过两层循环来处理数 i ：如果 i 未被标记为已访问（即 vis[i] == 0 ），将其认定为素数并存入 pm 数组，同时设置 mu[i] 为 -1 ， phi[i] 为 i - 1 。对于每个素数 pm[j] ，如果 i * pm[j] < N ，将 vis[i * pm[j]] 标记为已访问。如果 i 能被 pm[j] 整除，更新 phi[i * pm[j]] 并跳出内层循环；否则，更新 mu[i * pm[j]] 和 phi[i * pm[j]] 。最后计算 mu 数组和 phi 数组的前缀和。

函数 Sphi ：

• 首先进行快速判断，如果 n 小于预处理的范围 N ，直接返回预处理得到的 phi[n] 值。
• 如果 n 对应的结果未被计算过（即不在 mp_phi 映射中），则通过迭代计算：初始化 ans 为 n * (n + 1) / 2 。然后通过循环将 ans 减去子区间的结果乘以区间长度。
• 最后将计算结果存入 mp_phi 映射中并返回。

函数 Smu ：

• 与 Sphi 函数类似，先进行快速判断，如果 n 小于 N ，直接返回 mu[n] ；如果未计算过，通过迭代计算。
• 初始化 ans 为 1 ，通过循环更新 ans ，并将结果存入 mp_mu 映射中返回。

主函数 main ：调用 init 函数进行初始化。读入测试用例的数量 T 。在每个测试用例中，读入数 n ，分别调用 Sphi 和 Smu 函数计算并输出相应的结果。

P3768 简单的数学题

题目描述

P3768 简单的数学题 - 洛谷 | 计算机科学教育新生态 (luogu.com.cn)

运行代码

#include<bits/stdc++.h>
#define N 5000000
#define LL long long

using namespace std;

int vis[N], p[N], phi[N], cnt;
LL P, n, sum[N], inv;
map<LL, LL> mp;

// 计算 a 的 b 次幂对 P 取模
LL powerMod(LL a, LL b) {
    LL res = 1;
    while (b) {
        if (b & 1) res = res * a % P;
        a = a * a % P;
        b >>= 1;
    }
    return res;
}

// 初始化 phi 数组和 sum 数组
void init() {
    inv = powerMod(6, P - 2);
    phi[1] = 1;
    for (int i = 2; i < N; i++) {
        if (!vis[i]) {
            p[++cnt] = i;
            phi[i] = i - 1;
        }
        for (int j = 1; i * p[j] < N; j++) {
            vis[i * p[j]] = 1;
            if (i % p[j] == 0) {
                phi[i * p[j]] = p[j] * phi[i];
                break;
            }
            phi[i * p[j]] = (p[j] - 1) * phi[i];
        }
    }
    for (int i = 1; i < N; i++)
        sum[i] = (sum[i - 1] + 1LL * i * i % P * phi[i] % P) % P;
}

// 计算 x 的平方和
LL squareSum(LL x) {
    x %= P;
    return x * (x + 1) % P * (2 * x + 1) % P * inv % P;
}

// 计算 x 的立方和
LL cubeSum(LL x) {
    x %= P;
    return (x * (x + 1) / 2) % P * ((x * (x + 1) / 2) % P) % P;
}

// 计算 Sn
LL Sn(LL x) {
    if (x < N) return sum[x];
    if (mp[x]) return mp[x];
    LL res = cubeSum(x);
    for (LL l = 2, r; l <= x; l = r + 1) {
        r = x / (x / l);
        res = (res - (squareSum(r) - squareSum(l - 1)) % P * Sn(x / l) % P + P) % P;
    }
    return mp[x] = res;
}

// 计算最终结果
LL calculate() {
    LL ans = 0;
    for (LL l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l);
        ans = (ans + (Sn(r) - Sn(l - 1)) % P * cubeSum(n / l) % P + P) % P;
    }
    return ans;
}

int main() {
    scanf("%lld%lld", &P, &n);
    init();
    printf("%lld\n", calculate());
    return 0;
}

代码思路

1. 定义了一些常量和数组：N 用于定义一些计算的上限。vis 数组用于标记数字是否已被处理。p 数组用于存储质数。phi 数组用于存储欧拉函数值。sum 数组用于存储某些中间计算结果。

2. 定义了一些函数：qpow 函数用于快速计算幂取模运算。init 函数用于初始化 phi 数组和 sum 数组。在初始化 phi 数组时，使用了埃氏筛法的思想，同时计算出每个数的欧拉函数值。S2 函数用于计算平方和并取模。S3 函数用于计算立方和并取模。Sn 函数用于计算特定形式的和。如果 x < N ，直接返回预计算的 sum[x] ；否则通过分块计算的方法来求解。calc 函数用于最终的计算，通过迭代和调用前面定义的函数来得到结果。

3. 在 main 函数中，读取输入的 P 和 n ，调用 init 函数进行初始化，然后调用 calc 函数计算结果并输出。

P3803 【模板】多项式乘法（FFT）

题目描述

P3803 【模板】多项式乘法（FFT） - 洛谷 | 计算机科学教育新生态 (luogu.com.cn)

运行代码

#include <iostream>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
const int N = 4e6;
const double PI = acos(-1);

struct complex {
    double x, y;
    complex operator+(const complex& t) const {
        return {x + t.x, y + t.y};
    }
    complex operator-(const complex& t) const {
        return {x - t.x, y - t.y};
    }
    complex operator*(const complex& t) const {
        return {x * t.x - y * t.y, x * t.y + y * t.x};
    }
} A[N], B[N];
int R[N];

// 交换两个复数
void swapComplex(complex& a, complex& b) {
    complex temp = a;
    a = b;
    b = temp;
}

// 快速傅里叶变换
void FFT(complex A[], int n, int op) {
    for (int i = 0; i < n; ++i)
        R[i] = R[i / 2] / 2 + ((i & 1)? n / 2 : 0);

    for (int i = 0; i < n; ++i)
        if (i < R[i])
            swapComplex(A[i], A[R[i]]);

    for (int i = 2; i <= n; i <<= 1) {
        complex w1({cos(2 * PI / i), sin(2 * PI / i) * op});
        for (int j = 0; j < n; j += i) {
            complex wk({1, 0});
            for (int k = j; k < j + i / 2; ++k) {
                complex x = A[k], y = A[k + i / 2] * wk;
                A[k] = x + y;
                A[k + i / 2] = x - y;
                wk = wk * w1;
            }
        }
    }
}

int main() {
    int n, m;
    scanf("%d%d", &n, &m);
    for (int i = 0; i <= n; i++) scanf("%lf", &A[i].x);
    for (int i = 0; i <= m; i++) scanf("%lf", &B[i].x);
    for (m = n + m, n = 1; n <= m; n <<= 1);
    FFT(A, n, 1);
    FFT(B, n, 1);
    for (int i = 0; i < n; ++i) A[i] = A[i] * B[i];
    FFT(A, n, -1);
    for (int i = 0; i <= m; ++i) 
        printf("%d ", (int)(A[i].x / n + 0.5));
    return 0;
}

代码思路

FFT 是一种用于快速计算离散傅里叶变换（Discrete Fourier Transform，DFT）的算法。其核心思想是利用分治法，将长度为 n 的 DFT 分解为两个长度为 n/2 的 DFT，从而大大降低计算复杂度。

1. 首先定义了一个 complex 结构体来表示复数，其中包含实部 x 和虚部 y 。
2. FFT 函数用于执行快速傅里叶变换。第一步通过计算 R 数组来确定交换元素的位置。然后进行位逆序交换，将输入的复数序列按照特定的顺序重新排列。接下来通过循环，对于每个长度为 i 的子区间进行合并计算。其中 w1 是旋转因子，wk 用于在每个子区间内逐步乘以旋转因子进行计算。
3. 在 main 函数中，首先读取两个多项式的系数。然后通过扩展长度进行 FFT 计算，将两个多项式在频域相乘，再通过逆 FFT 转换回时域，得到最终的乘积结果。