library

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

---

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

Mobius table
(Math/mobius.hpp)

• View this file on GitHub
• Last update: 2024-12-26 05:48:15+09:00
• Include: #include "Math/mobius.hpp"

Depends on

• Divisor Multiple Transform (Convolution/divisor.hpp)
• Prime Sieve (Math/sieve.hpp)

Code

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

vector<int> Mobius(ll n) {
    vector<int> mu(n + 1);
    mu[1] = 1;
    DivisorTransform::mobius(mu);
    return mu;
}

/**
 * @brief Mobius table
 */

#line 2 "Math/sieve.hpp"

template<int L=50101010>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/mobius.hpp"

vector<int> Mobius(ll n) {
    vector<int> mu(n + 1);
    mu[1] = 1;
    DivisorTransform::mobius(mu);
    return mu;
}

/**
 * @brief Mobius table
 */

Back to top page