library

This documentation is automatically generated by online-judge-tools/verification-helper


Project maintained by tko919 Hosted on GitHub Pages — Theme by mattgraham

:warning: lpf table
(Math/totient.hpp)

Depends on

Code

#pragma once
#include "Convolution/divisor.hpp"

vector<int> phitable(int n) {
    vector<int> ret(n + 1);
    iota(ALL(ret), 0);
    DivisorTransform::mobius(ret);
    return ret;
}

/**
 * @brief lpf table
 */
#line 2 "Math/sieve.hpp"

template <int L = 101010101> vector<int> sieve(int N) {
    bitset<L> isp;
    int n, sq = ceil(sqrt(N));
    for (int z = 1; z <= 5; z += 4) {
        for (int y = z; y <= sq; y += 6) {
            for (int x = 1; x <= sq and (n = 4 * x * x + y * y) <= N; ++x) {
                isp[n].flip();
            }
            for (int x = y + 1; x <= sq and (n = 3 * x * x - y * y) <= N;
                 x += 2) {
                isp[n].flip();
            }
        }
    }
    for (int z = 2; z <= 4; z += 2) {
        for (int y = z; y <= sq; y += 6) {
            for (int x = 1; x <= sq and (n = 3 * x * x + y * y) <= N; x += 2) {
                isp[n].flip();
            }
            for (int x = y + 1; x <= sq and (n = 3 * x * x - y * y) <= N;
                 x += 2) {
                isp[n].flip();
            }
        }
    }
    for (int y = 3; y <= sq; y += 6) {
        for (int z = 1; z <= 2; ++z) {
            for (int x = z; x <= sq and (n = 4 * x * x + y * y) <= N; x += 3) {
                isp[n].flip();
            }
        }
    }
    for (int n = 5; n <= sq; ++n)
        if (isp[n]) {
            for (int k = n * n; k <= N; k += n * n) {
                isp[k] = false;
            }
        }
    isp[2] = isp[3] = true;

    vector<int> ret;
    for (int i = 2; i <= N; i++)
        if (isp[i]) {
            ret.push_back(i);
        }
    return ret;
}

/**
 * @brief Prime Sieve
 */
#line 3 "Convolution/divisor.hpp"

namespace DivisorTransform{
    int n;
    vector<int> ps;
    template<typename T>void zeta(vector<T>& a){
        int N=a.size()-1;
        if(n<N){
            ps=sieve(N);
            n=N;
        }
        for(auto& p:ps){
            for(int k=1;k*p<=N;k++)a[k*p]+=a[k];
        }
    }
    template<typename T>void mobius(vector<T>& a){
        int N=a.size()-1;
        if(n<N){
            ps=sieve(N);
            n=N;
        }
        for(auto& p:ps){
            for(int k=N/p;k;k--)a[k*p]-=a[k];
        }
    }
};

namespace MultipleTransform{
    int n;
    vector<int> ps;
    template<typename T>void zeta(vector<T>& a){
        int N=a.size()-1;
        if(n<N){
            ps=sieve(N);
            n=N;
        }
        for(auto& p:ps){
            for(int k=N/p;k;k--)a[k]+=a[k*p];
        }
    }
    template<typename T>void mobius(vector<T>& a){
        int N=a.size()-1;
        if(n<N){
            ps=sieve(N);
            n=N;
        }
        for(auto& p:ps){
            for(int k=1;k*p<=N;k++)a[k]-=a[k*p];
        }
    }
};

/**
 * @brief Divisor Multiple Transform
 * @docs docs/divisor.md
 */
#line 3 "Math/totient.hpp"

vector<int> phitable(int n) {
    vector<int> ret(n + 1);
    iota(ALL(ret), 0);
    DivisorTransform::mobius(ret);
    return ret;
}

/**
 * @brief lpf table
 */
Back to top page