結果

問題 No.2620 Sieve of Coins
ユーザー 👑 NachiaNachia
提出日時 2024-01-15 01:01:19
言語 C++17(gcc12)
(gcc 12.3.0 + boost 1.87.0)
結果
AC  
実行時間 73 ms / 2,000 ms
コード長 6,024 bytes
コンパイル時間 1,536 ms
コンパイル使用メモリ 114,964 KB
実行使用メモリ 11,136 KB
最終ジャッジ日時 2024-09-28 02:09:03
合計ジャッジ時間 4,303 ms
ジャッジサーバーID
（参考情報） 		judge1 / judge5
このコードへのチャレンジ
（要ログイン）
ファイルパターン 結果
sample AC * 5
other AC * 53
権限があれば一括ダウンロードができます

ソースコード

diff # 

// 00:32:43

#ifdef NACHIA
#define _GLIBCXX_DEBUG
#else
#define NDEBUG
#endif
#include <iostream>
#include <string>
#include <vector>
#include <algorithm>
#include <utility>
#include <queue>
#include <array>
#include <cmath>


#include <cassert>


namespace nachia{

namespace prime_sieve_explicit_internal{
    std::vector<bool> isprime = { false }; // a[x] := isprime(2x+1)

    void CalcIsPrime(int z){
        if((int)isprime.size() *2+1 < z+1){
            int new_z = isprime.size();
            while(new_z*2+1 < z+1) new_z *= 2;
            z = new_z-1;
            isprime.resize(z+1, true);
            for(int i=1; i*(i+1)*2<=z; i++) if(isprime[i]){
                for(int j=i*(i+1)*2; j<=z; j+=i*2+1) isprime[j] = false;
            }
        }
    }
    
    std::vector<int> prime_list = {2};
    int prime_list_max = 0;

    void CalcPrimeList(int z){
        while((int)prime_list.size() < z){
            if((int)isprime.size() <= prime_list_max + 1) CalcIsPrime(prime_list_max * 2 + 10);
            for(int p=prime_list_max+1; p<(int)isprime.size(); p++){
                if(isprime[p]) prime_list.push_back(p*2+1);
            }
            prime_list_max = isprime.size() - 1;
        }
    }

    void CalcPrimeListUntil(int z){
        if(prime_list_max < z){
            CalcIsPrime(z);
            for(int p=prime_list_max+1; p<(int)isprime.size(); p++){
                if(isprime[p]) prime_list.push_back(p*2+1);
            }
            prime_list_max = isprime.size() - 1;
        }
    }

}


bool IsprimeExplicit(int n){
    using namespace prime_sieve_explicit_internal;
    if(n == 2) return true;
    if(n % 2 == 0) return false;
    CalcIsPrime(n);
    return isprime[(n-1)/2];
}

int NthPrimeExplicit(int n){
    using namespace prime_sieve_explicit_internal;
    CalcPrimeList(n);
    return prime_list[n];
}

int PrimeCountingExplicit(int n){
    using namespace prime_sieve_explicit_internal;
    if(n < 2) return 0;
    CalcPrimeListUntil(n);
    auto res = std::upper_bound(prime_list.begin(), prime_list.end(), n) - prime_list.begin();
    return (int)res;
}

// [l, r)
std::vector<bool> SegmentedSieveExplicit(long long l, long long r){
    assert(0 <= l); assert(l <= r);
    long long d = r - l;
    if(d == 0) return {};
    std::vector<bool> res(d, true);
    for(long long p=2; p*p<=r; p++) if(IsprimeExplicit(p)){
        long long il = (l+p-1)/p, ir = (r+p-1)/p;
        if(il <= p) il = p;
        for(long long i=il; i<ir; i++) res[i*p-l] = false;
    }
    if(l < 2) for(long long p=l; p<2 && p<r; p++) res[l-p] = false;
    return res;
}


} // namespace nachia


namespace nachia{

template<class Elem>
void DivisorZeta(std::vector<Elem>& a){
    using namespace prime_sieve_explicit_internal;
    int n = a.size() - 1;
    for(int d=2; d<=n; d++) if(IsprimeExplicit(d)) for(int i=1; i*d<=n; i++) a[i*d] += a[i];
}

template<class Elem>
void DivisorInvZeta(std::vector<Elem>& a){
    using namespace prime_sieve_explicit_internal;
    int n = a.size() - 1;
    for(int d=2; d<=n; d++) if(IsprimeExplicit(d)) for(int i=n/d; i>=1; i--) a[i] += a[i*d];
}

template<class Elem>
void DivisorMobius(std::vector<Elem>& a){
    using namespace prime_sieve_explicit_internal;
    int n = a.size() - 1;
    for(int d=2; d<=n; d++) if(IsprimeExplicit(d)) for(int i=n/d; i>=1; i--) a[i*d] -= a[i];
}

template<class Elem>
void DivisorInvMobius(std::vector<Elem>& a){
    using namespace prime_sieve_explicit_internal;
    int n = a.size() - 1;
    for(int d=2; d<=n; d++) if(IsprimeExplicit(d)) for(int i=1; i*d<=n; i++) a[i] -= a[i*d];
}

template<class Elem>
std::vector<Elem> GcdConvolution(std::vector<Elem> a, std::vector<Elem> b){
    assert(a.size() == b.size());
    assert(1 <= a.size());
    DivisorInvZeta(a);
    DivisorInvZeta(b);
    for(int i=1; i<(int)a.size(); i++) a[i] *= b[i];
    DivisorInvMobius(a);
    return a;
}

template<class Elem>
std::vector<Elem> LcmConvolution(std::vector<Elem> a, std::vector<Elem> b){
    assert(a.size() == b.size());
    assert(1 <= a.size());
    DivisorZeta(a);
    DivisorZeta(b);
    for(int i=1; i<(int)a.size(); i++) a[i] *= b[i];
    DivisorMobius(a);
    return a;
}

}
using namespace std;
using i64 = long long;
using u64 = unsigned long long;
#define rep(i,n) for(i64 i=0; i<(i64)(n); i++)
#define repr(i,n) for(i64 i=(i64)(n)-1; i>=0; i--)
const i64 INF = 1001001001001001001;
const char* yn(bool x){ return x ? "Yes" : "No"; }
template<typename A> void chmin(A& l, const A& r){ if(r < l) l = r; }
template<typename A> void chmax(A& l, const A& r){ if(l < r) l = r; }
template<typename A> using nega_queue = priority_queue<A,vector<A>,greater<A>>;

int main(){
    i64 L, N; cin >> L >> N;
    vector<i64> A(N); rep(i,N) cin >> A[i];

    vector<vector<i64>> pow23(1);
    pow23[0].push_back(1);
    while(pow23.back()[0] * 3 <= L){
        i64 l = pow23.back()[0];
        pow23.emplace_back(1);
        pow23.back()[0] = l * 3;
    }
    for(auto& a : pow23) while(a.back() * 2 <= L) a.push_back(a.back() * 2);
    int H = (int)max(pow23.size(), pow23[0].size());
    vector<vector<int>> coins(H, vector<int>(H));
    for(i64 a : A){
        int p2 = 0; while(a%2 == 0){ p2++; a /= 2; }
        int p3 = 0; while(a%3 == 0){ p3++; a /= 3; }
        coins[p3][p2] = 1;
    }
    rep(y,H) rep(x,H) if(coins[y][x]){
        for(int yy=y; yy<H; yy++) for(int xx=x; xx<H; xx++) if(xx != x || yy != y) coins[yy][xx] ^= 1;
    }
    i64 ans = 0;
    rep(y,H) rep(x,H) if(coins[y][x]) if(y<(int)pow23.size()) if(x<(int)pow23[y].size()){
        i64 mp = pow23[y][x];
        i64 smallL = L / mp;
        i64 sz = 1; while((sz+1)*(sz+1) <= smallL) sz++;
        vector<i64> Q(sz+1);
        for(i64 i=1; i<=sz; i++) if(i%2 != 0 && i%3 != 0){
            Q[i] += smallL / (i * i);
            Q[i] -= smallL / (i * i) / 2;
            Q[i] -= smallL / (i * i) / 3;
            Q[i] += smallL / (i * i) / 6;
        }
        nachia::DivisorInvMobius(Q);
        ans += Q[1];
    }
    cout << ans << endl;

    return 0;
}
0