結果

問題 No.2266 Fractions (hard)
ユーザー 👑 NachiaNachia
提出日時 2023-04-08 14:56:42
言語 C++17
(gcc 12.3.0 + boost 1.83.0)
結果
AC  
実行時間 84 ms / 6,000 ms
コード長 11,889 bytes
コンパイル時間 838 ms
コンパイル使用メモリ 89,584 KB
実行使用メモリ 5,760 KB
最終ジャッジ日時 2024-10-03 10:04:27
合計ジャッジ時間 2,695 ms
ジャッジサーバーID
（参考情報） 		judge3 / judge1
このコードへのチャレンジ
（要ログイン）

テストケース

テストケース表示
入力 結果 実行時間
実行使用メモリ
testcase_00 AC 3 ms
5,248 KB
testcase_01 AC 60 ms
5,632 KB
testcase_02 AC 34 ms
5,760 KB
testcase_03 AC 2 ms
5,248 KB
testcase_04 AC 2 ms
5,248 KB
testcase_05 AC 2 ms
5,248 KB
testcase_06 AC 6 ms
5,248 KB
testcase_07 AC 4 ms
5,248 KB
testcase_08 AC 2 ms
5,248 KB
testcase_09 AC 1 ms
5,248 KB
testcase_10 AC 2 ms
5,248 KB
testcase_11 AC 1 ms
5,248 KB
testcase_12 AC 2 ms
5,248 KB
testcase_13 AC 59 ms
5,376 KB
testcase_14 AC 6 ms
5,248 KB
testcase_15 AC 7 ms
5,248 KB
testcase_16 AC 54 ms
5,248 KB
testcase_17 AC 44 ms
5,248 KB
testcase_18 AC 62 ms
5,248 KB
testcase_19 AC 56 ms
5,504 KB
testcase_20 AC 6 ms
5,248 KB
testcase_21 AC 69 ms
5,632 KB
testcase_22 AC 33 ms
5,248 KB
testcase_23 AC 7 ms
5,760 KB
testcase_24 AC 7 ms
5,504 KB
testcase_25 AC 72 ms
5,632 KB
testcase_26 AC 8 ms
5,632 KB
testcase_27 AC 59 ms
5,632 KB
testcase_28 AC 8 ms
5,760 KB
testcase_29 AC 7 ms
5,632 KB
testcase_30 AC 8 ms
5,632 KB
testcase_31 AC 7 ms
5,632 KB
testcase_32 AC 59 ms
5,632 KB
testcase_33 AC 2 ms
5,248 KB
testcase_34 AC 38 ms
5,504 KB
testcase_35 AC 39 ms
5,632 KB
testcase_36 AC 84 ms
5,632 KB
testcase_37 AC 76 ms
5,632 KB
権限があれば一括ダウンロードができます

ソースコード

diff # 

#line 2 "D:\\Programming\\VSCode\\competitive-cpp\\nachia\\math\\floor-sum.hpp"

#include <cassert>
#include <utility>

namespace nachia{


// a : any value
// mod != 0
std::pair<long long, unsigned long long> SafeDiv(long long a, unsigned long long mod){
    using u64 = unsigned long long;
    if(a >= 0) return std::make_pair(0, (u64)a);
    if(mod >= (u64)1 << 62) return std::make_pair(-1, (u64)a + mod);
    long long q = a / mod;
    long long m = a % (long long)mod;
    if(m){ q--; m += mod; }
    return std::make_pair(q, m);
}

unsigned long long nC2Uint64(unsigned long long n){
    return (n%2) ? ((n-1)/2*n) : (n/2*(n-1));
}

// n : any
// 1 <= m
// a : any
// b : any
// n * a%m + b%m < 2**64
unsigned long long FloorSumU64Unsigned(
    unsigned long long n,
    unsigned long long m,
    unsigned long long a,
    unsigned long long b
){
    using u64 = unsigned long long;
    assert(1 <= m);
    u64 ans = 0;
    while(n){
        if(a >= m){
            ans += a / m * nC2Uint64(n);
            a %= m;
        }
        if(b >= m){
            ans += b / m * n;
            b %= m;
        }
        u64 y_max = a * n + b;
        if (y_max < m) return ans;
        n = y_max / m;
        b = y_max % m;
        y_max = a; a = m; m = y_max;
    }
    return ans;
}

// n : any
// 1 <= m
// a : any
// b : any
// (n+1) * m < 2**64
unsigned long long FloorSumU64Signed(
    unsigned long long n,
    unsigned long long m,
    long long a,
    long long b
){

    using u64 = unsigned long long;
    auto ua = SafeDiv(a, m);
    auto ub = SafeDiv(b, m);
    u64 ans = FloorSumU64Unsigned(n, m, ua.second, ub.second);
    ans += ua.first / m * nC2Uint64(n);
    ans += ub.first / m * n;
    return ans;
}

} // namespace nachia
#line 4 "D:\\Programming\\VSCode\\competitive-cpp\\nachia\\math\\rational-number-search.hpp"

namespace nachia{

class RationalNumberSearch{
public:
    RationalNumberSearch(unsigned long long maxVal){
        assert(maxVal < 1ull << 63);
        mx = maxVal;
    }
    bool hasNext(){ return state >= 0; }
    std::pair<unsigned long long, unsigned long long> getNext() const {
        switch(state){
        case 0: return { a0+a1, b0+b1 };
        case 1: return { a0+tr*a1, b0+tr*b1 };
        case 2: return { a1+tr*a0, b1+tr*b0 };
        case 3: return { a0+(tl+tr)/2*a1, b0+(tl+tr)/2*b1 };
        case 4: return { a1+(tl+tr)/2*a0, b1+(tl+tr)/2*b0 };
        }
        return {0,0};
    }
    void give(bool toRight){
        int x = toRight ? 1 : 0;
        switch(state){
        case 0:
            tl = 1; tr = 2;
            if(a0 + a1 > mx || b0 + b1 > mx){ state = -1; }
            else{ state = (toRight ? 1 : 2); }
            break;
        case 1: case 2:
            if(x ^ (2-state)){ state += 2; }
            else{ tr *= 2; tl *= 2; }
            break;
        case 3: case 4:
            ((x ^ (4-state)) ? tr : tl) = (tl+tr)/2;
            break;
        }
        while(givecheck());
    }
private:
    using UInt = unsigned long long;
    UInt a0=0, b0=1, a1=1, b1=0, tl=0, tr=0, mx;
    int state = 0;
    bool givecheck(){
        auto st = [this](int x){ state = x; return true; };
        auto trq = [this](UInt x0, UInt x1) -> bool {
            bool f = x0+tr*x1 > mx;
            if(f) tr = (mx-x0)/x1 + 1;
            return f;
        };
        bool f = false;
        switch(state){
        case -1 : break;
        case 0:
            if(a0 + a1 > mx || b0 + b1 > mx){ state = -1; }
            break;
        case 1:
            if(trq(a0,a1)) f = true;
            if(trq(b0,b1)) f = true;
            if(f) return st(3);
            break;
        case 2:
            if(trq(a1,a0)) f = true;
            if(trq(b1,b0)) f = true;
            if(f) return st(4);
            break;
        case 3:
            if(tl + 1 == tr){
                a0 += a1 * tl;
                b0 += b1 * tl;
                return st(0);
            }
            break;
        case 4:
            if(tl + 1 == tr){
                a1 += a0 * tl;
                b1 += b0 * tl;
                return st(0);
            }
            break;
        }
        return false;
    }
};

} // namespace nachia
#line 2 "D:\\Programming\\VSCode\\competitive-cpp\\nachia\\math\\multiplicative-convolution-usefloat.hpp"

#line 4 "D:\\Programming\\VSCode\\competitive-cpp\\nachia\\math\\multiplicative-convolution-usefloat.hpp"
#include <vector>
#include <algorithm>

namespace nachia{

// https://maspypy.com/dirichlet-%E7%A9%8D%E3%81%A8%E3%80%81%E6%95%B0%E8%AB%96%E9%96%A2%E6%95%B0%E3%81%AE%E7%B4%AF%E7%A9%8D%E5%92%8C#toc10
template<class Val>
struct MultiplicativeConvolution{
    using MyType = MultiplicativeConvolution;
    static long long Div(long long a, long long d){
        return (long long)((double)a / (double)d);
    }
    struct Status{
        long long N, K, L;
        Status(long long N){
            this->N = N;
            if(N <= 10){ K=N; L=1; }
            else if(N <= 1000){
                K = 1; while(K*K < N) K++;
                L = (N+K-1) / K;
            }
            else{
                L = 1; while(L*L*L/50 < N) L++;
                K = (N+L-1) / L;
            }
        }
        bool operator==(Status r) const {
            return N == r.N && K == r.K && L == r.L;
        }
    };
    Status stat;
    std::vector<Val> a;
    std::vector<Val> A;
    MultiplicativeConvolution(Status give_stat) : stat(give_stat), a(stat.K+1, 0), A(stat.L+1, 0){}
    static MyType Zeta(Status stat){
        MyType res(stat);
        for(long long i=1; i<=stat.K; i++) res.a[i] = Val(1);
        for(long long i=1; i<=stat.L; i++) res.A[i] = Val(Div(stat.N,i));
        return res;
    }
    static MyType Zeta1(Status stat){
        MyType res(stat);
        Val inv2 = Val(1) / Val(2);
        for(long long i=1; i<=stat.K; i++) res.a[i] = Val(i);
        for(long long i=1; i<=stat.L; i++) res.A[i] = Val(Div(stat.N,i)) * Val(Div(stat.N,i)+1) * inv2;
        return res;
    }
    static MyType One(Status stat){
        MyType res(stat);
        res.a[1] = Val(1);
        for(long long i=1; i<=stat.L; i++) res.A[i] = Val(1);
        return res;
    }
    MyType& operator+=(const MyType& r) {
        assert(this->stat == r.stat);
        for(long long i=1; i<=stat.K; i++) a[i] += r.a[i];
        for(long long i=1; i<=stat.L; i++) A[i] += r.A[i];
        return *this;
    }
    MyType operator+(const MyType& r) const { auto res = *this; return res += r; }
    MyType& operator-=(const MyType& r) {
        assert(this->stat == r.stat);
        for(long long i=1; i<=stat.K; i++) a[i] -= r.a[i];
        for(long long i=1; i<=stat.L; i++) A[i] -= r.A[i];
        return *this;
    }
    MyType operator-(const MyType& r) const { auto res = *this; return res -= r; }
    MyType operator*(const MyType& r) const { return mult(*this,r); }
    MyType& operator*=(const MyType& r) { return *this = mult(*this,r); }
    MyType operator/(const MyType& r) const { return div(*this,r); }
    MyType& operator/=(const MyType& r) { return *this = div(*this,r); }
private:
    static MyType mult(const MyType& l, const MyType& r){
        assert(l.stat == r.stat);
        long long N = l.stat.N, K = l.stat.K, L = l.stat.L;
        std::vector<Val> lFairy(K+1, 0);
        for(long long i=1; i<=K; i++) lFairy[i] = lFairy[i-1] + l.a[i];
        std::vector<Val> rFairy(K+1, 0);
        for(long long i=1; i<=K; i++) rFairy[i] = rFairy[i-1] + r.a[i];
        MyType res(l.stat);
        for(long long i=1; i<=K; i++) for(long long j=1; i*j<=K; j++) res.a[i*j] += l.a[i] * r.a[j];
        for(long long c=L,C=Div(N,L),M=1; c>=1; --c){ // calc res.A[c]
            C = Div(N,c);
            while(M*M<C) M++; // M = ceil sqrt N/c
            Val AM = (M<=K ? lFairy[M] : l.A[Div(N,M)]);
            for(long long i=1; i<=M; i++) res.A[c] += l.a[i] * (c*i<=L ? r.A[c*i] : rFairy[Div(C,i)]);
            for(long long j=1; M*j<=C; j++) res.A[c] += r.a[j] * ((c*j<=L ? l.A[c*j] : lFairy[Div(C,j)]) - AM);
        }
        return res;
    }
    static MyType div(const MyType& l, const MyType& r){
        Val AInv = Val(1) / Val(r.a[1]);
        assert(l.stat == r.stat);
        long long N = l.stat.N, K = l.stat.K, L = l.stat.L;
        MyType res(l.stat);
        for(long long i=1; i<=K; i++){
            res.a[i] = (l.a[i] - res.a[i]) * AInv;
            for(long long j=2; i*j<=K; j++) res.a[i*j] += res.a[i] * r.a[j];
        }
        std::vector<Val> resFairy(K+1, 0);
        for(long long i=1; i<=K; i++) resFairy[i] = resFairy[i-1] + res.a[i];
        std::vector<Val> rFairy(K+1, 0);
        for(long long i=1; i<=K; i++) rFairy[i] = rFairy[i-1] + r.a[i];
        for(long long c=L,C=Div(N,L),M=1; c>=1; --c){ // calc res.A[c]
            C = Div(N,c);
            while(M*M<C) M++; // M = ceil sqrt N/c
            Val AM = rFairy[M];
            for(long long i=2; i<=M; i++) res.A[c] += r.a[i] * (c*i<=L ? res.A[c*i] : resFairy[Div(C,i)]);
            for(long long j=1; M*j<=C; j++) res.A[c] += res.a[j] * ((c*j<=L ? r.A[c*j] : rFairy[Div(C,j)]) - AM);
            res.A[c] = (l.A[c] - res.A[c]) * AInv;
        }
        return res;
    }
};

} // namespace nachia

#line 4 "D:\\Programming\\VSCode\\competitive-cpp\\nachia\\math\\enumerate-quotients.hpp"

namespace nachia{

template<class Int = unsigned long long>
std::vector<Int> EnumerateQuotients(Int N, bool doInsertZero){
    std::vector<Int> res;
    if(doInsertZero) res.push_back(0);
    for(Int k=1; k*k<N; k++) res.push_back(k);
    int qp1 = res.size();
    for(Int k=1; k*k<=N; k++) res.push_back(N/k);
    std::reverse(res.begin() + qp1, res.end());
    return res;
}

}
#line 5 "..\\Main.cpp"
#include <iostream>
#include <string>
#line 9 "..\\Main.cpp"
using namespace std;
using i64 = long long;
#define rep(i,n) for(int i=0; i<(int)(n); i++)
const i64 INF = 1001001001001001001;

pair<i64,i64> testcase(){
    i64 N, K; cin >> N >> K;

    vector<i64> quotients = nachia::EnumerateQuotients<i64>(N, true);
    int numQ = quotients.size()-1;

    using Mf = nachia::MultiplicativeConvolution<i64>;
    auto mfstat = Mf::Status(N);
    Mf zeta = Mf::Zeta(mfstat);
    Mf mobius = Mf::One(mfstat) / zeta;
    vector<i64> mertensRaw(mfstat.K+1);
    for(i64 i=1; i<=mfstat.K; i++) mertensRaw[i] = mertensRaw[i-1] + mobius.a[i];

    auto getMertens = [&](i64 n) -> i64 {
        if(n <= mfstat.K) return mertensRaw[n];
        return mobius.A[N/n];
    };

    auto countFracs = [&]() -> i64 {
        i64 fracCnt = 0;
        rep(q,numQ){
            i64 times = getMertens(quotients[q+1]) - getMertens(quotients[q]);
			i64 v = quotients[numQ-q];
            i64 cnt = v * (v+1) / 2;
            fracCnt += times * cnt;
        }
        return fracCnt;
    };
    
    i64 cntall = countFracs();
    if(cntall * 2 - 1 < K) return {-1,-1};
    if(cntall == K) return {1,1};
    bool sw = false;
    if(cntall < K){ K = cntall*2 - K; sw = true; }

    auto srch = nachia::RationalNumberSearch(N);
    i64 a = 0, b = 0;
	i64 sqrtN = 1; while(sqrtN*sqrtN <= N) sqrtN++;
	sqrtN++;
	sqrtN *= 4;
    while(srch.hasNext()){
        i64 ax = srch.getNext().first;
        i64 bx = srch.getNext().second;
        i64 fracCnt = 0;
		vector<i64> cntn(sqrtN+1);
		for(i64 i=1; i<=sqrtN; i++) cntn[i] = cntn[i-1] + i * ax / bx;
        rep(q,numQ){
            i64 times = getMertens(quotients[q+1]) - getMertens(quotients[q]);
			i64 v = quotients[numQ-q];
            i64 cnt = v <= sqrtN ? cntn[v] : nachia::FloorSumU64Signed(quotients[numQ-q]+1, bx, ax, 0);
            fracCnt += times * cnt;
        }
        if(K <= fracCnt){ a = ax; b = bx; }
        srch.give(fracCnt < K);
    }

    if(sw) swap(a, b);
    return {a,b};
}

int main(){
	ios::sync_with_stdio(false);
	cin.tie(nullptr);
    auto ans = testcase();
    if(ans.first < 0) cout << "-1\n"; else cout << ans.first << '/' << ans.second << '\n';
    return 0;
}
0