GCD/LCM Convolution
(convolution/gcd_lcm_convolution.hpp)

• View this file on GitHub
• Last update: 2024-01-07 22:37:45+09:00
• Include: #include "convolution/gcd_lcm_convolution.hpp"

Depends on

• Multiple/Divisor Fast Zeta/Möbius Transform (convolution/divisor_zeta_moebius_transform.hpp)

Verified with

• test/yosupo/gcd_convolution.test.cpp
• test/yosupo/lcm_convolution.test.cpp

Code

#pragma once
#include <vector>

#include "divisor_zeta_moebius_transform.hpp"

/**
 * @brief GCD/LCM Convolution
 */

template <typename T>
std::vector<T> gcd_convolution(std::vector<T> a, std::vector<T> b) {
    const int n = std::max(a.size(), b.size());
    a.resize(n);
    b.resize(n);
    multiple_fzt(a);
    multiple_fzt(b);
    for (int i = 0; i < n; ++i) a[i] *= b[i];
    multiple_fmt(a);
    return a;
}

template <typename T>
std::vector<T> lcm_convolution(std::vector<T> a, std::vector<T> b) {
    const int n = std::max(a.size(), b.size());
    a.resize(n);
    b.resize(n);
    divisor_fzt(a);
    divisor_fzt(b);
    for (int i = 0; i < n; ++i) a[i] *= b[i];
    divisor_fmt(a);
    return a;
}

#line 2 "convolution/gcd_lcm_convolution.hpp"
#include <vector>

#line 3 "convolution/divisor_zeta_moebius_transform.hpp"

template <typename T>
void divisor_fzt(std::vector<T>& a) {
    const int n = a.size();
    std::vector<bool> sieve(n, true);
    for (int p = 2; p < n; ++p) {
        if (!sieve[p]) continue;
        for (int k = 1; k * p < n; ++k) {
            sieve[k * p] = false;
            a[k * p] += a[k];
        }
    }
}

template <typename T>
void divisor_fmt(std::vector<T>& a) {
    const int n = a.size();
    std::vector<bool> sieve(n, true);
    for (int p = 2; p < n; ++p) {
        if (!sieve[p]) continue;
        for (int k = (n - 1) / p; k > 0; --k) {
            sieve[k * p] = false;
            a[k * p] -= a[k];
        }
    }
}

template <typename T>
void multiple_fzt(std::vector<T>& a) {
    const int n = a.size();
    std::vector<bool> sieve(n, true);
    for (int p = 2; p < n; ++p) {
        if (!sieve[p]) continue;
        for (int k = (n - 1) / p; k > 0; --k) {
            sieve[k * p] = false;
            a[k] += a[k * p];
        }
    }
}

template <typename T>
void multiple_fmt(std::vector<T>& a) {
    const int n = a.size();
    std::vector<bool> sieve(n, true);
    for (int p = 2; p < n; ++p) {
        if (!sieve[p]) continue;
        for (int k = 1; k * p < n; ++k) {
            sieve[k * p] = false;
            a[k] -= a[k * p];
        }
    }
}
#line 5 "convolution/gcd_lcm_convolution.hpp"

/**
 * @brief GCD/LCM Convolution
 */

template <typename T>
std::vector<T> gcd_convolution(std::vector<T> a, std::vector<T> b) {
    const int n = std::max(a.size(), b.size());
    a.resize(n);
    b.resize(n);
    multiple_fzt(a);
    multiple_fzt(b);
    for (int i = 0; i < n; ++i) a[i] *= b[i];
    multiple_fmt(a);
    return a;
}

template <typename T>
std::vector<T> lcm_convolution(std::vector<T> a, std::vector<T> b) {
    const int n = std::max(a.size(), b.size());
    a.resize(n);
    b.resize(n);
    divisor_fzt(a);
    divisor_fzt(b);
    for (int i = 0; i < n; ++i) a[i] *= b[i];
    divisor_fmt(a);
    return a;
}

Back to top page