結果
| 問題 |
No.2266 Fractions (hard)
|
| ユーザー |
 Nachia
|
| 提出日時 | 2023-04-08 14:56:42 |
| 言語 | C++17 (gcc 13.3.0 + boost 1.87.0) |
| 結果 |
AC
|
| 実行時間 | 86 ms / 6,000 ms |
| コード長 | 11,889 bytes |
| コンパイル時間 | 880 ms |
| コンパイル使用メモリ | 86,804 KB |
| 最終ジャッジ日時 | 2025-02-12 03:53:28 |
|
ジャッジサーバーID （参考情報） |
judge1 / judge1 |
（要ログイン）
| ファイルパターン | 結果 |
|---|---|
| sample | AC * 3 |
| other | AC * 35 |
ソースコード
#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;
}