結果
| 問題 |
No.2718 Best Consonance
|
| コンテスト | |
| ユーザー |
tonegawa
|
| 提出日時 | 2024-04-05 22:33:36 |
| 言語 | C++17 (gcc 13.3.0 + boost 1.87.0) |
| 結果 |
TLE
|
| 実行時間 | - |
| コード長 | 49,975 bytes |
| コンパイル時間 | 3,267 ms |
| コンパイル使用メモリ | 188,940 KB |
| 最終ジャッジ日時 | 2025-02-20 21:46:47 |
|
ジャッジサーバーID (参考情報) |
judge2 / judge3 |
(要ログイン)
| ファイルパターン | 結果 |
|---|---|
| sample | AC * 4 |
| other | AC * 10 TLE * 1 -- * 25 |
ソースコード
#include <iostream>
#include <string>
#include <vector>
#include <array>
#include <tuple>
#include <stack>
#include <queue>
#include <deque>
#include <algorithm>
#include <set>
#include <map>
#include <unordered_set>
#include <unordered_map>
#include <bitset>
#include <cmath>
#include <functional>
#include <cassert>
#include <climits>
#include <iomanip>
#include <numeric>
#include <memory>
#include <random>
#include <thread>
#include <chrono>
#define allof(obj) (obj).begin(), (obj).end()
#define range(i, l, r) for(int i=l;i<r;i++)
#define unique_elem(obj) obj.erase(std::unique(allof(obj)), obj.end())
#define bit_subset(i, S) for(int i=S, zero_cnt=0;(zero_cnt+=i==S)<2;i=(i-1)&S)
#define bit_kpop(i, n, k) for(int i=(1<<k)-1,x_bit,y_bit;i<(1<<n);x_bit=(i&-i),y_bit=i+x_bit,i=(!i?(1<<n):((i&~y_bit)/x_bit>>1)|y_bit))
#define bit_kth(i, k) ((i >> k)&1)
#define bit_highest(i) (i?63-__builtin_clzll(i):-1)
#define bit_lowest(i) (i?__builtin_ctzll(i):-1)
#define sleepms(t) std::this_thread::sleep_for(std::chrono::milliseconds(t))
using ll = long long;
using ld = long double;
using ul = uint64_t;
using pi = std::pair<int, int>;
using pl = std::pair<ll, ll>;
using namespace std;
template<typename F, typename S>
std::ostream &operator<<(std::ostream &dest, const std::pair<F, S> &p){
dest << p.first << ' ' << p.second;
return dest;
}
template<typename T>
std::ostream &operator<<(std::ostream &dest, const std::vector<std::vector<T>> &v){
int sz = v.size();
if(sz==0) return dest;
for(int i=0;i<sz;i++){
int m = v[i].size();
for(int j=0;j<m;j++) dest << v[i][j] << (i!=sz-1&&j==m-1?'\n':' ');
}
return dest;
}
template<typename T>
std::ostream &operator<<(std::ostream &dest, const std::vector<T> &v){
int sz = v.size();
if(sz==0) return dest;
for(int i=0;i<sz-1;i++) dest << v[i] << ' ';
dest << v[sz-1];
return dest;
}
template<typename T, size_t sz>
std::ostream &operator<<(std::ostream &dest, const std::array<T, sz> &v){
if(sz==0) return dest;
for(int i=0;i<sz-1;i++) dest << v[i] << ' ';
dest << v[sz-1];
return dest;
}
template<typename T>
std::ostream &operator<<(std::ostream &dest, const std::set<T> &v){
for(auto itr=v.begin();itr!=v.end();){
dest << *itr;
itr++;
if(itr!=v.end()) dest << ' ';
}
return dest;
}
template<typename T, typename E>
std::ostream &operator<<(std::ostream &dest, const std::map<T, E> &v){
for(auto itr=v.begin();itr!=v.end();){
dest << '(' << itr->first << ", " << itr->second << ')';
itr++;
if(itr!=v.end()) dest << '\n';
}
return dest;
}
std::ostream &operator<<(std::ostream &dest, __int128_t value) {
std::ostream::sentry s(dest);
if (s) {
__uint128_t tmp = value < 0 ? -value : value;
char buffer[128];
char *d = std::end(buffer);
do {
--d;
*d = "0123456789"[tmp % 10];
tmp /= 10;
} while (tmp != 0);
if (value < 0) {
--d;
*d = '-';
}
int len = std::end(buffer) - d;
if (dest.rdbuf()->sputn(d, len) != len) {
dest.setstate(std::ios_base::badbit);
}
}
return dest;
}
template<typename T>
vector<T> make_vec(size_t sz, T val){return std::vector<T>(sz, val);}
template<typename T, typename... Tail>
auto make_vec(size_t sz, Tail ...tail){
return std::vector<decltype(make_vec<T>(tail...))>(sz, make_vec<T>(tail...));
}
template<typename T>
vector<T> read_vec(size_t sz){
std::vector<T> v(sz);
for(int i=0;i<(int)sz;i++) std::cin >> v[i];
return v;
}
template<typename T, typename... Tail>
auto read_vec(size_t sz, Tail ...tail){
auto v = std::vector<decltype(read_vec<T>(tail...))>(sz);
for(int i=0;i<(int)sz;i++) v[i] = read_vec<T>(tail...);
return v;
}
long long max(long long a, int b){return std::max(a, (long long)b);}
long long max(int a, long long b){return std::max((long long)a, b);}
long long min(long long a, int b){return std::min(a, (long long)b);}
long long min(int a, long long b){return std::min((long long)a, b);}
long long modulo(long long a, long long m){a %= m; return a < 0 ? a + m : a;}
template<typename T>
struct safe_vector : std::vector<T>{
using std::vector<T>::vector;
T& operator [](size_t i){return this->at(i);}
};
template<typename T, int N>
struct safe_array : std::array<T, N>{
using std::array<T, N>::array;
T& operator [](size_t i){return this->at(i);}
};
ll ceil_div(ll x, ll y){
assert(y > 0);
return (x + (x > 0 ? y - 1 : 0)) / y;
}
ll floor_div(ll x, ll y){
assert(y > 0);
return (x + (x > 0 ? 0 : -y + 1)) / y;
}
void io_init(){
std::cin.tie(nullptr);
std::ios::sync_with_stdio(false);
}
#include <cstdint>
struct barrett_reduction{
unsigned int mod;
unsigned long long m;
barrett_reduction(unsigned int _mod) : mod(_mod){
m = ((__uint128_t)1 << 64) / mod;
}
unsigned int reduce(unsigned int a){
unsigned long long q = ((__uint128_t)a * m) >> 64;
a -= q * mod; // 0 <= a < 2 * mod
// return a;
return a >= mod ? a - mod : a;
}
unsigned int mul(unsigned int a, unsigned int b){
return reduce((unsigned long long)a * b);
}
// {gcd(a, mod), x}, such that a * x ≡ gcd(a, mod)
std::pair<unsigned int, unsigned int> inv(unsigned int a){
if(a >= mod) a = reduce(a);
if(a == 0) return {mod, 0};
unsigned int s = mod, t = a;
long long m0 = 0, m1 = 1;
while(t){
int u = s / t;
s -= t * u;
m0 -= m1 * u;
std::swap(m0, m1);
std::swap(s, t);
}
if(m0 < 0) m0 += mod / s;
return {s, m0};
}
};
// 64bit mod対応
// R = 2^64
// 偶数modだと壊れる
struct montgomery_reduction_64bit{
private:
// [0, 2 * MOD)
inline uint64_t reduce_unsafe(__uint128_t x) const{
x = (x + ((uint64_t)x * (uint64_t)NEG_INV) * MOD) >> 64;
return x;
}
void _set_mod(uint64_t mod){
assert((mod < (1ULL << 63)) && (mod & 1));
MOD = mod;
NEG_INV = 0;
__uint128_t s = 1, t = 0;
for(int i = 0; i < 64; i++){
if (~t & 1) {
t += MOD;
NEG_INV += s;
}
t >>= 1;
s <<= 1;
}
R2 = ((__uint128_t)1 << 64) % MOD;
R2 = R2 * R2 % MOD;
ONE = generate(1);
}
__uint128_t MOD, NEG_INV, R2;
uint64_t ONE;
public:
montgomery_reduction_64bit(){}
montgomery_reduction_64bit(uint64_t mod){_set_mod(mod);}
void set_mod(uint64_t mod){
_set_mod(mod);
}
uint64_t mod()const{
return MOD;
}
inline uint64_t one()const{
return ONE;
}
inline uint64_t generate(uint64_t x)const{
return reduce((__uint128_t)x * R2);
}
inline uint64_t reduce(__uint128_t x)const{
x = (x + ((uint64_t)x * (uint64_t)NEG_INV) * MOD) >> 64;
return x < MOD ? x : x - MOD;
}
inline uint64_t fix(uint64_t x)const{
return x < MOD ? x : x - MOD;
}
// [0, 2MOD)
inline uint64_t mul(uint64_t mx, uint64_t my)const{
return reduce_unsafe((__uint128_t)mx * my);
}
inline uint64_t mul_safe(uint64_t mx, uint64_t my)const{
return reduce((__uint128_t)mx * my);
}
// [0, 2MOD)
inline uint64_t add(uint64_t mx, uint64_t my)const{
return (mx >= MOD ? mx - MOD : mx) + (my >= MOD ? my - MOD : my);
}
inline uint64_t add_safe(uint64_t mx, uint64_t my)const{
uint64_t res = (mx >= MOD ? mx - MOD : mx) + (my >= MOD ? my - MOD : my);
return res >= MOD ? res - MOD : res;
}
// [0, 2MOD)
inline uint64_t sub(uint64_t mx, uint64_t my)const{
if(my >= MOD) my -= MOD;
return mx >= my ? mx - my : mx + MOD - my;
}
inline uint64_t sub_safe(uint64_t mx, uint64_t my)const{
if(my >= MOD) my -= MOD;
uint64_t res = mx >= my ? mx - my : mx + MOD - my;
return res >= MOD ? res - MOD : res;
}
inline uint64_t pow(uint64_t ma, uint64_t b)const{
uint64_t ret = one();
while(b){
if(b & 1) ret = mul(ret, ma);
ma = mul(ma, ma);
b >>= 1;
}
return ret;
}
inline uint64_t pow_safe(uint64_t ma, uint64_t b)const{
return fix(pow(ma, b));
}
};
unsigned long long mod_pow_mr(unsigned long long a, unsigned long long b, unsigned long long m){
montgomery_reduction_64bit mr(m);
return mr.reduce(mr.pow(mr.generate(a), b));
}
namespace prime_sieve{
std::vector<int> primes, min_factor;// 素数, 各数を割り切る最小の素数
// O(MAX_N loglog MAX_N)
// [1, MAX_N]を扱えるように初期化
void init(int MAX_N){
min_factor.resize(MAX_N + 1, -1);
for(int i = 2; i <= MAX_N; i++){
if(min_factor[i] == -1){
primes.push_back(i);
min_factor[i] = i;
}
for(int p : primes){
if((long long)p * i > MAX_N || p > min_factor[i]) break;
min_factor[p * i] = p;
}
}
}
bool is_prime(int n){
assert(n < min_factor.size());
return n == min_factor[n];
}
// {{素因数, 数}}, O(log n)
std::vector<std::pair<int, int>> factorize(int n){
assert(n < min_factor.size());
std::vector<std::pair<int, int>> ret;
while(n > 1){
int cnt = 0, f = min_factor[n];
while(n % f == 0){
n /= f;
cnt++;
}
ret.push_back({f, cnt});
}
return ret;
}
// 約数列挙, O(√n)
std::vector<int> divisor(int n){
auto p = factorize(n);
std::vector<std::vector<int>> x;
for(int i = 0; i < p.size(); i++){
x.push_back(std::vector<int>{1});
for(int j = 0; j < p[i].second; j++) x[i].push_back(x[i][j] * p[i].first);
}
int l = 0, r = 1;
std::vector<int> ret{1};
for(int i = 0; i < x.size(); i++){
for(auto e : x[i]){
for(int j = l; j < r; j++){
ret.push_back(ret[j] * e);
}
}
l = r;
r = ret.size();
}
return std::vector<int>(ret.begin() + l, ret.end());
}
// O(logN)
unsigned long long totient_function(unsigned long long n){
unsigned long long res = n;
int prev = -1;
while(n > 1){
if(min_factor[n] > prev){
res -= res / min_factor[n];
prev = min_factor[n];
}
n /= min_factor[n];
}
return res;
}
int mobius_function(int x){
int pcnt = 0;
while(x > 1){
int y = x / min_factor[x];
if(min_factor[x] == min_factor[y]) return 0;
x = y;
pcnt++;
}
return pcnt % 2 == 0 ? 1 : -1;
}
};
bool _miller_rabin_mr(unsigned long long n, const montgomery_reduction_64bit &mr){
static std::vector<int> small_p{2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47};
static std::vector<unsigned long long> A{2, 325, 9375, 28178, 450775, 9780504, 1795265022};
static std::vector<unsigned long long> B{2, 7, 61};
if(n <= 1) return false;
if(n <= 50){
for(int l = n < 20 ? 0 : 8, r = n < 20 ? 8 : 15; l < r; l++) if(small_p[l] == n) return true;
return false;
}
if(!(n & 1)) return false;
unsigned long long d = n - 1;
unsigned long long one = mr.one(), mone = mr.generate(n - 1);
d >>= __builtin_ctzll(d);
for(unsigned long long a : (n >> 32) ? A : B){
if(a % n == 0) continue;
unsigned long long d2s = d; // d * 2^s, 0 <= s <= (n - 1)が2で割れる回数
unsigned long long y = mr.pow_safe(mr.generate(a), d);
while(d2s != n - 1 && y != one && y != mone){
y = mr.mul_safe(y, y);
d2s <<= 1;
}
if(y != mone && !(d2s & 1)) return false;
}
return true;
}
bool miller_rabin_mr(unsigned long long n){
if(n % 2 == 0) return n == 2 ? true : false;
montgomery_reduction_64bit mr(n);
return _miller_rabin_mr(n, mr);
}
// https://en.wikipedia.org/wiki/Binary_GCD_algorithm
unsigned long long binary_gcd(unsigned long long a, unsigned long long b){
if(!a || !b) return !a ? b : a;
int shift = __builtin_ctzll(a | b); // [1] gcd(2a', 2b') = 2 * gcd(a', b')
a >>= __builtin_ctzll(a);
do{
// if b is odd
// gcd(2a', b) = gcd(a', b), if a = 2a'(a is even)
// gcd(a, b) = gcd(|a - b|, min(a, b)), if a is odd
b >>= __builtin_ctzll(b); // make b odd
if(a > b) std::swap(a, b);
b -= a;
}while(b); // gcd(a, 0) = a
return a << shift; // [1]
}
unsigned long long generate_random_prime(unsigned long long min_n = 2, unsigned long long max_n = ~0ULL){
std::random_device seed_gen;
std::mt19937_64 engine(seed_gen());
__uint128_t len = max_n - min_n + 1;
// https://en.wikipedia.org/wiki/Prime_number_theorem
while(true){
unsigned long long a = engine() % len + min_n;
if(miller_rabin_mr(a)){
return a;
}
}
}
namespace rho_factorization{
unsigned long long rho(unsigned long long n){
static std::vector<int> small_p{2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47};
for(int sp : small_p) if(n % sp == 0) return sp; // n < 50
montgomery_reduction_64bit mr(n);
if(_miller_rabin_mr(n, mr)) return n;
auto try_factorize = [n, mr](unsigned long long c){
c = mr.generate(c);
auto f = [mr, c](unsigned long long mx){
return mr.add(mr.mul(mx, mx), c);
};
unsigned long long m = 1LL << ((64 - __builtin_clzll(n)) / 8);
unsigned long long y = n, r = 1, q = 1;
unsigned long long x, g, k, ys;
do{
x = y;
y = mr.generate(y);
for(int i = 0; i < r; i++) y = f(y);
y = mr.reduce(y);
k = 0;
while(k < r && g == 1){
q = mr.generate(q);
y = mr.generate(y);
ys = y;
for(int i = 0; i < std::min(m, r - k); i++){
y = f(y);
unsigned long long z = mr.reduce(y);
q = mr.mul(q, mr.generate(x > z ? x - z : z - x));
}
y = mr.reduce(y);
g = binary_gcd(mr.reduce(q), n);
k += m;
}
r <<= 1;
}while(g == 1);
if(g == n){
do{
ys = f(ys);
unsigned long long z = mr.reduce(ys);
g = binary_gcd(x > z ? x - z : z - x, n);
}while(g == 1);
}
return g; // g == n when failure
};
unsigned long long c = 1, res = n;
do{
res = try_factorize(c);
// c = generate_random_prime(2, n - 1);
c = (c + 1) % n;
}while(res == n);
return res;
}
std::vector<unsigned long long> factorize(unsigned long long n){
if(n <= 1) return {};
unsigned long long x = rho(n);
if(x == n) return {x};
auto l = factorize(x);
auto r = factorize(n / x);
l.insert(l.end(), r.begin(), r.end());
return l;
}
// {素数, 個数}の形で返す
std::vector<std::pair<unsigned long long, int>> factorize2(unsigned long long n){
auto p = factorize(n);
sort(p.begin(), p.end());
std::vector<std::pair<unsigned long long, int>> ret;
for(int i : p){
if(ret.empty() || ret.back().first != i) ret.push_back({i, 1});
else ret.back().second++;
}
return ret;
}
// 素因数の集合(重複なし, ソート済)を返す
std::vector<unsigned long long> prime_factor(unsigned long long n){
auto p = factorize(n);
std::sort(p.begin(), p.end());
p.erase(std::unique(p.begin(), p.end()), p.end());
return p;
}
// 10^18以下の高度合成数 897612484786617600の約数が103680個なので全列挙して良さそう
std::vector<unsigned long long> divisor(unsigned long long n){
auto p = factorize(n);
std::sort(p.begin(), p.end());
std::vector<std::pair<unsigned long long, int>> x;
for(int i = 0; i < p.size(); i++){
if(!i || p[i] != p[i - 1]) x.push_back({p[i], 1});
else x.back().second++;
}
int sz = 1;
for(auto [p_cur, cnt] : x) sz *= cnt + 1;
std::vector<unsigned long long> res(sz);
res[0] = 1;
int r_prev = 1, r = 1;
for(auto [p_cur, cnt] : x){
unsigned long long ppow = 1;
for(int c = 0; c < cnt; c++){
ppow *= p_cur;
for(int i = 0; i < r_prev; i++){
res[r++] = res[i] * ppow;
}
}
r_prev = r;
}
return res;
}
int mobius_function(long long x){
auto P = rho_factorization::factorize(x);
for(long long p : P) if((x / p) % p == 0) return 0;
return P.size() % 2 == 0 ? 1 : -1;
}
unsigned long long totient_function(unsigned long long n){
unsigned long long res = n;
auto prims = rho_factorization::prime_factor(n);
for(auto p : prims) res -= res / p;
return res;
}
};
// p: 素数
unsigned long long find_primitive_root(unsigned long long p){
static std::random_device seed_gen;
static std::mt19937_64 engine(seed_gen());
//assert(miller_rabin_mr(p));
auto primes = rho_factorization::prime_factor(p - 1);
while(true){
bool f = true;
unsigned long long a = engine() % (p - 1) + 1;
for(unsigned long long pk : primes){
// a ^ (p - 1) / pk ≡ 1 (mod p) -> no
if(mod_pow_mr(a, (p - 1) / pk, p) == 1){
f = false;
break;
}
}
if(f) return a;
}
}
struct bigint_prime_factor{
int sq, max_x;
std::vector<int> low;
std::unordered_map<int, int> high;
bigint_prime_factor(int max_x): sq(sqrt(max_x)), max_x(max_x), low(sq + 1, 0){
// 篩が初期化済か
assert(prime_sieve::min_factor.size() > max_x);
}
// O(log(x))
// x^kを掛ける
void mul(int x, int k = 1){
assert(0 < x && x <= max_x);
while(x > 1){
int p = prime_sieve::min_factor[x];
if(p <= sq) low[p] += k;
else{
auto [itr, f] = high.emplace(x, k);
if(!f) itr->second += k;
}
x /= p;
}
}
// O(log(x))
// x^kで割る(割り切れない場合エラー)
void div(int x, int k = 1){
assert(0 < x && x <= max_x);
while(x > 1){
int p = prime_sieve::min_factor[x];
if(p <= sq){
assert(low[p] >= k);
low[p] -= k;
}else{
auto itr = high.find(p);
assert(itr != high.end() && itr->second >= k);
itr->second -= k;
if(itr->second == 0) high.erase(itr);
}
x /= p;
}
}
// 素数pで何回割れるか, O(1)
// pが素数でない場合0
int k_divisible_prime(int p){
if(p > max_x || !prime_sieve::is_prime(p)) return 0;
if(p > sq){
auto itr = high.find(p);
return itr == high.end() ? 0 : itr->second;
}
return low[p];
}
// xで何回割れるか, O(log(x))
// x = 1のときinf回割れる
// @param 0 < x <= 篩の最大値
int k_divisible(int x){
assert(x > 0);
static constexpr int inf = std::numeric_limits<int>::max();
int ans = inf;
for(auto [p, num] : prime_sieve::factorize(x)){
if(p > sq){
auto itr = high.find(p);
if(itr == high.end()) return 0;
ans = std::min(ans, itr->second / num);
}else{
ans = std::min(ans, low[p] / num);
}
}
return ans;
}
// rで何回割り切れるか, O(sqrt(r.max_x)以下の素数 + rの素因数の種類)
int k_divisible(const bigint_prime_factor &r){
static constexpr int inf = std::numeric_limits<int>::max();
int ans = inf;
for(int i = 0; prime_sieve::primes[i] <= r.sq; i++){
int p = prime_sieve::primes[i];
if(!r.low[p]) continue;
if(p <= sq){
if(low[p] < r.low[p]) return 0;
ans = std::min(ans, low[p] / r.low[p]);
}else{
auto itr = high.find(p);
if(itr->second < r.low[p]) return 0;
ans = std::min(ans, itr->second / r.low[p]);
}
}
for(auto [p, num] : r.high){
assert(num);
if(p <= sq){
if(low[p] < num) return 0;
ans = std::min(ans, low[p] / num);
}else{
auto itr = high.find(p);
if(itr->second < num) return 0;
ans = std::min(ans, itr->second / num);
}
}
return ans;
}
};
// res[i] := v[i]の位数
template<typename mint>
std::vector<long long> multiplicative_order_many(const std::vector<mint> &v){
int n = v.size();
std::vector<long long> res(n, 1);
for(auto [p, q] : rho_factorization::factorize2(mint::mod() - 1)){
long long x = mint(p).pow(q).val(), y = (mint::mod() - 1) / x;
for(int i = 0; i < n; i++){
long long z = x;
for(int j = 0; j < q; j++){
z /= p;
if(v[i].pow(y * z).val() != 1){
res[i] *= z * p;
break;
}
}
}
}
return res;
}
template<typename Val>
struct binary_indexed_tree{
int M, H;
std::vector<Val> sum;
binary_indexed_tree(){}
binary_indexed_tree(int N): M(N), H(31 - __builtin_clz(M)), sum(M + 1 , 0){}
binary_indexed_tree(const std::vector<Val> &v): M(v.size()), H(31 - __builtin_clz(M)), sum(1){
sum.insert(sum.begin() + 1, v.begin(), v.end());
for(int i = 1; i <= M; i++){
int nxt = i + (i & (-i));
if(nxt <= M) sum[nxt] += sum[i];
}
}
void update(int k, Val x){
for(int i = k + 1; i <= M; i += (i & (-i))) sum[i] += x;
}
Val query(int r){
Val ret = 0;
for(int k = r; k > 0; k -= (k & (-k))) ret += sum[k];
return ret;
}
Val query(int l, int r){
return query(r) - query(l);
}
// sum[0, k]がx以上になるような最小のkとsum[0, k], 無い場合は{M, 総和}
// sumが単調非減少であることが必要
using p = std::pair<int, Val>;
p lower_bound(Val x){
int v = 1 << H, h = H;
Val s = 0, t = 0;
while(h--){
if(M < v) v -= 1 << h;
else if(x <= s + sum[v]) t = s + sum[v], v -= 1 << h;
else s += sum[v], v += 1 << h;
}
if(v == M + 1) return {M, s};
return (x <= s + sum[v] ? p{v - 1, s + sum[v]} : p{v, t});
}
};
template<typename Idx, typename Val>
struct binary_indexed_tree_compressed{
int M;
std::vector<std::pair<Val, Val>> sum;
binary_indexed_tree_compressed(){}
binary_indexed_tree_compressed(int N): M(N), sum(N + 1, {0, 0}){}
binary_indexed_tree_compressed(const std::vector<Val> &v): M(v.size()), sum(M + 1, {0, 0}){
for(int i = 0; i < M; i++) sum[i + 1].first = v[i];
for(int i = 1; i <= M; i++){
int nxt = i + (i & (-i));
if(nxt <= M) sum[nxt].first += sum[i].first;
}
}
void update(int kc, Val x){
for(int i = kc + 1; i <= M; i += (i & (-i))) sum[i].first += x;
}
// [lc, rc) (元の数直線上では[l, r))にzを加算
void update(Idx l, Idx r, int lc, int rc, Val z){
Val a = l * z, b = r * z;
for(int i = lc + 1; i <= M; i += (i & (-i))){
sum[i].first -= a;
sum[i].second += z;
}
for(int i = rc + 1; i <= M; i += (i & (-i))){
sum[i].first += b;
sum[i].second -= z;
}
}
Val query(Idx r, int rc){
Val a = 0, b = 0;
for(int i = rc; i > 0; i -= (i & (-i))){
a += sum[i].first;
b += sum[i].second;
}
return a + (b * r);
}
// [lc, rc) (元の数直線上では[l, r))の和
Val query(Idx l, Idx r, int lc, int rc){
if(lc >= rc) return 0;
return query(r, rc) - query(l, lc);
}
};
#include <limits>
template<typename Key, typename Val>
struct map_sum_cache{
static constexpr Key inf = std::numeric_limits<Key>::max();
static constexpr Val inf_val = std::numeric_limits<Val>::max();
static constexpr int limit_size_per_node = 16;
private:
struct node{
int h, sz, sz_sum;
std::array<Key, limit_size_per_node> keys;
std::array<Val, limit_size_per_node> vals;
Val sum, sumsub;
node *l, *r;
node(): h(1), l(nullptr), r(nullptr), sum(0){}
node(Key _key, Val _val): h(1), sz(1), sz_sum(1), l(nullptr), r(nullptr){keys[0] = _key; vals[0] = sum = sumsub = _val;}
node(const std::vector<std::pair<Key, Val>> &v, int l, int r): h(1), sz(r - l), sz_sum(sz), sum(0), l(nullptr), r(nullptr){
assert(sz < limit_size_per_node);
for(int i = 0; i < sz; i++){
keys[i] = v[l + i].first;
vals[i] = v[l + i].second;
sum += vals[i];
}
sumsub = sum;
}
int balanace_factor(){
return (l ? l->h : 0) - (r ? r->h : 0);
}
node *split_half(){
assert(sz == limit_size_per_node);
node *u = new node();
sz = limit_size_per_node / 2;
u->sz_sum = u->sz = limit_size_per_node - sz;
for(int i = 0; i < u->sz; i++){
u->keys[i] = keys[sz + i];
u->vals[i] = vals[sz + i];
u->sum += u->vals[i];
}
u->sumsub = u->sum;
sum -= u->sum;
return u;
}
};
node *root;
int size(node *v){return v ? v->sz_sum : 0;}
void update(node *v){
v->h = std::max(v->l ? v->l->h : 0, v->r ? v->r->h : 0) + 1;
v->sz_sum = (v->l ? v->l->sz_sum : 0) + (v->r ? v->r->sz_sum : 0) + v->sz;
v->sumsub = (v->l ? v->l->sumsub : 0) + (v->r ? v->r->sumsub : 0) + v->sum;
}
node *rotate_right(node *v){
node *l = v->l;
v->l = l->r;
l->r = v;
update(v);
update(l);
return l;
}
node *rotate_left(node *v){
node *r = v->r;
v->r = r->l;
r->l = v;
update(v);
update(r);
return r;
}
node *balance(node *v){
int bf = v->balanace_factor();
assert(-2 <= bf && bf <= 2);
if(bf == 2){
if(v->l->balanace_factor() == -1){
v->l = rotate_left(v->l);
update(v);
}
return rotate_right(v);
}else if(bf == -2){
if(v->r->balanace_factor() == 1){
v->r = rotate_right(v->r);
update(v);
}
return rotate_left(v);
}
return v;
}
node *build(const std::vector<node*> &nodes, int l, int r){
int m = (l + r) >> 1;
node *v = nodes[m];
if(m > l) v->l = build(nodes, l, m);
if(r > m + 1) v->r = build(nodes, m + 1, r);
update(v);
return v;
}
node *insert_leftmost(node *v, node *u){
if(!v) return u;
v->l = insert_leftmost(v->l, u);
update(v);
return balance(v);
}
node *emplace_inner(node *v, Key k, Val val, bool replace = false){
if(!v) return new node(k, val);
if(v->l && k < v->keys[0]){
v->l = emplace_inner(v->l, k, val, replace);
}else if(v->r && v->keys[v->sz - 1] < k){
v->r = emplace_inner(v->r, k, val, replace);
}else{
int idx = std::lower_bound(v->keys.begin(), v->keys.begin() + v->sz, k) - v->keys.begin();
if(idx < v->sz && v->keys[idx] == k){
v->vals[idx] += val;
v->sum += val;
v->sumsub += val;
return v;
}
for(int i = v->sz; i > idx; i--){
v->keys[i] = v->keys[i - 1];
v->vals[i] = v->vals[i - 1];
}
v->keys[idx] = k;
v->vals[idx] = val;
v->sum += val;
v->sz++;
if(v->sz == limit_size_per_node){
v->r = insert_leftmost(v->r, v->split_half());
}
}
update(v);
return balance(v);
}
Val query_inner(node *v, Key k){
Val ret = 0;
while(v){
if(k < v->keys[0]){
v = v->l;
}else if(k > v->keys[v->sz - 1]){
ret += (v->l ? v->l->sumsub : 0) + v->sum;
v = v->r;
}else{
ret += (v->l ? v->l->sumsub : 0);
for(int i = 0; i < v->sz; i++){
if(v->keys[i] >= k) return ret;
ret += v->vals[i];
}
}
}
return ret;
}
public:
map_sum_cache(): root(nullptr){}
map_sum_cache(std::vector<std::pair<Key, Val>> v){
std::sort(v.begin(), v.end());
init_sorted(v);
}
// すでにソート済みかつキーがユニーク
void init_sorted(const std::vector<std::pair<Key, Val>> &v){
if(v.empty()){
root = nullptr;
return;
}
int n = v.size();
int size_per_node = limit_size_per_node / 2;
int m = (n + size_per_node - 1) / size_per_node;
std::vector<node*> nodes(m);
for(int i = 0; i < m; i++){
nodes[i] = new node(v, i * size_per_node, std::min((i + 1) * size_per_node, n));
}
root = build(nodes, 0, m);
}
int size(){
return size(root);
}
void update(Key key, Val val){
root = emplace_inner(root, key, val);
}
Val query(Key l, Key r){
return query_inner(root, r) - query_inner(root, l);
}
};
template<typename Idx, typename Val>
struct offline_static_rectangle_sum{
private:
struct Point{
Idx x, y;
Val z;
};
struct Query{
Idx lx, rx, ly, ry;
};
std::vector<Point> P;
std::vector<Query> Q;
public:
void update(Idx x, Idx y, Val z){P.push_back(Point{x, y, z});}
void query(Idx lx, Idx rx, Idx ly, Idx ry){Q.push_back(Query{lx, rx, ly, ry});}
std::vector<Val> solve(){
struct Event{
Idx x;
int ly, ry;
int id;
};
int N = Q.size();
if(P.empty() || Q.empty()) return std::vector<Val>(N, 0);
std::vector<Event> Q2(2 * N);
std::vector<Val> ans(N, 0);
std::vector<Idx> Y;
std::sort(P.begin(), P.end(), [](const Point &a, const Point &b){return a.y < b.y;});
for(Point &t : P){
if(Y.empty() || Y.back() != t.y) Y.push_back(t.y);
t.y = int(Y.size()) - 1;
}
for(int i = 0; i < N; i++){
int ly = std::lower_bound(Y.begin(), Y.end(), Q[i].ly) - Y.begin();
int ry = std::lower_bound(Y.begin(), Y.end(), Q[i].ry) - Y.begin();
Q2[2 * i] = Event{Q[i].lx, ly, ry, i};
Q2[2 * i + 1] = Event{Q[i].rx, ly, ry, i + N};
}
std::sort(P.begin(), P.end(), [](const Point &a, const Point &b){return a.x < b.x;});
std::sort(Q2.begin(), Q2.end(), [](const Event &a, const Event &b){return a.x < b.x;});
binary_indexed_tree<Val> bit(Y.size());
int p = 0, q = 0;
while(q < 2 * N){
if(p == P.size() || Q2[q].x <= P[p].x){
Val sum = (!p ? Val(0) : bit.query(Q2[q].ly, Q2[q].ry));
if(Q2[q].id >= N) ans[Q2[q].id - N] += sum;
else ans[Q2[q].id] -= sum;
q++;
}else{
bit.update(P[p].y, P[p].z);
p++;
}
}
return ans;
}
};
template<typename Idx, typename Val>
struct offline_point_add_rectangle_sum{
private:
struct query_type{
int id;
Idx lx, rx, ly, ry;
Val z;
query_type(int a, Idx b, Idx c, Idx d, Idx e, Val f): id(a), lx(b), rx(c), ly(d), ry(e), z(f){}
};
struct event_type{
Idx x;
int q, id;
int lyc, ryc;
Val z;
event_type(Idx a, int b, int c, int d, int e, Val f): x(a), q(b), id(c), lyc(d), ryc(e), z(f){}
};
std::vector<query_type> q;
std::vector<int> qcount{0};
int qs = 0;
void solve(int l, int r, std::vector<Val> &ans){
if(r - l < 2) return;
int mid = (l + r) >> 1;
solve(l, mid, ans);
solve(mid, r, ans);
int left_update = (mid - l) - (qcount[mid] - qcount[l]);
int right_query= qcount[r] - qcount[mid];
if(left_update == 0 || right_query == 0) return;
// compress y
std::vector<Idx> y;
for(int i = l; i < mid; i++) if(q[i].id == -1) y.push_back(q[i].ly);
std::sort(y.begin(), y.end());
y.erase(std::unique(y.begin(), y.end()), y.end());
binary_indexed_tree<Val> bit(y.size());
std::vector<event_type> e;
for(int i = l; i < mid; i++){
if(q[i].id == -1){
int y_idx = std::lower_bound(y.begin(), y.end(), q[i].ly) - y.begin();
e.push_back(event_type(q[i].lx, 2, -1, y_idx, 0, q[i].z));
}
}
for(int i = mid; i < r; i++){
if(q[i].id != -1){
int ly_idx = std::lower_bound(y.begin(), y.end(), q[i].ly) - y.begin();
int ry_idx = std::lower_bound(y.begin(), y.end(), q[i].ry) - y.begin();
e.push_back(event_type(q[i].lx, 0, q[i].id, ly_idx, ry_idx, 0));
e.push_back(event_type(q[i].rx, 1, q[i].id, ly_idx, ry_idx, 0));
}
}
std::sort(e.begin(), e.end(), [](event_type &a, event_type &b){
if(a.x != b.x) return a.x < b.x;
return a.q < b.q;
});
bool updated = false;
for(event_type &ei : e){
if(ei.q == 0){
if(updated) ans[ei.id] -= bit.query(ei.lyc, ei.ryc);
}else if(ei.q == 1){
if(updated) ans[ei.id] += bit.query(ei.lyc, ei.ryc);
}else{
updated = true;
bit.update(ei.lyc, ei.z);
}
}
}
public:
offline_point_add_rectangle_sum(){}
void update(Idx x, Idx y, Val z){
q.push_back(query_type(-1, x, 0, y, 0, z));
qcount.push_back(0);
}
void query(Idx lx, Idx rx, Idx ly, Idx ry){
q.push_back(query_type(qs++, lx, rx, ly, ry, 0));
qcount.push_back(1);
}
std::vector<Val> solve(){
std::vector<Val> ret(qs, 0);
for(int i = 1; i < qcount.size(); i++) qcount[i] += qcount[i - 1];
solve(0, q.size(), ret);
return ret;
}
};
template<typename Idx, typename Val>
struct partial_offline_point_add_rectangle_sum{
private:
bool is_built = false;
int M = 1;
std::vector<Idx> X;
std::vector<std::vector<Idx>> Y;
std::vector<binary_indexed_tree<Val>> BITs;
using point = std::tuple<Idx, Idx, Val>;
int lbx(Idx x){
return std::lower_bound(X.begin(), X.end(), x) - X.begin();
}
public:
std::vector<point> P;
void set_initial_point(Idx x, Idx y, Val z){
assert(!is_built);
P.push_back({x, y, z});
}
void build(){
assert(!is_built);
is_built = true;
int n = P.size();
std::sort(P.begin(), P.end(), [](const point &a, const point &b){
return std::get<1>(a) < std::get<1>(b);
});
for(int i = 0; i < n; i++) X.push_back(std::get<0>(P[i]));
std::sort(X.begin(), X.end());
X.erase(std::unique(X.begin(), X.end()), X.end());
M = X.size();
std::vector<std::vector<Val>> tmp(M + 1);
BITs.resize(M + 1);
Y.resize(M + 1);
for(int i = 0; i < n; i++){
auto [x, y, z] = P[i];
int k = lbx(x);
for(int j = k + 1; j <= M; j += (j & (-j))){
if(Y[j].empty()||Y[j].back()!=y){
Y[j].push_back(y);
tmp[j].push_back(z);
}else{
tmp[j].back() += z;
}
}
}
for(int i = 0; i <= M; i++) BITs[i] = binary_indexed_tree<Val>(tmp[i]);
}
partial_offline_point_add_rectangle_sum(){}
partial_offline_point_add_rectangle_sum(const std::vector<point> &v){
P = v;
build();
}
void update(Idx x, Idx y, Val z){
assert(is_built);
int xidx = lbx(x);
for(int i = xidx + 1; i <= M; i += (i & (-i))){
auto yidx = std::lower_bound(Y[i].begin(), Y[i].end(), y) - Y[i].begin();
BITs[i].update(yidx, z);
}
}
Val query(int rx, Idx ly, Idx ry){
assert(is_built);
Val ret = 0;
for(int i = rx; i > 0; i -= (i & (-i))){
int ridx = std::lower_bound(Y[i].begin(), Y[i].end(), ry) - Y[i].begin();
int lidx = std::lower_bound(Y[i].begin(), Y[i].end(), ly) - Y[i].begin();
ret += BITs[i].query(lidx, ridx);
}
return ret;
}
Val query(Idx lx, Idx rx, Idx ly, Idx ry){
return query(lbx(rx), ly, ry) - query(lbx(lx), ly, ry);
}
};
template<typename Idx, typename Val, int max_n_log = 30>
struct online_point_add_rectangle_sum{
private:
struct node{
map_sum_cache<Idx, Val> mp;
std::vector<std::pair<Idx, Val>> initial_points; // 初期化用
node *l, *r;
node(): l(nullptr), r(nullptr){}
};
node *root;
void update_inner(Idx x, Idx y, Val z){
node *v = root;
Idx lx = 0, rx = (Idx)1 << max_n_log;
while(true){
Idx mx = (lx + rx) / 2;
if(rx <= x){
if(!v->r) v->r = new node();
v = v->r;
lx = rx, rx += rx - mx;
}else{
v->mp.update(y, z);
if(rx - 1 == x) return;
rx = mx;
if(!v->l) v->l = new node();
v = v->l;
}
}
}
Val query_inner(Idx x, Idx ly, Idx ry){
Idx lx = 0, rx = (Idx)1 << max_n_log;
Val ret = 0;
node *v = root;
while(v){
Idx mx = (lx + rx) / 2;
if(rx <= x){
ret += v->mp.query(ly, ry);
if(rx == x) return ret;
v = v->r;
lx = rx;
rx += rx - mx;
}else{
v = v->l;
rx = mx;
}
}
return ret;
}
using point = std::tuple<Idx, Idx, Val>;
public:
online_point_add_rectangle_sum(): root(new node()){}
online_point_add_rectangle_sum(std::vector<point> v): root(new node()){
sort(v.begin(), v.end(), [](const point &a, const point &b){
return std::get<1>(a) < std::get<1>(b);
});
auto push = [&](Idx x, Idx y, Val z){
node *v = root;
Idx lx = 0, rx = (Idx)1 << max_n_log;
while(true){
Idx mx = (lx + rx) / 2;
if(rx <= x){
if(!v->r) v->r = new node();
v = v->r;
lx = rx, rx += rx - mx;
}else{
if(v->initial_points.empty() || v->initial_points.back().first != y){
v->initial_points.push_back({y, z});
}else{
v->initial_points.back().second += z;
}
if(rx - 1 == x) return;
rx = mx;
if(!v->l) v->l = new node();
v = v->l;
}
}
};
for(auto [x, y, z] : v) push(x, y, z);
auto init = [&](auto &&init, node *v) -> void {
v->mp.init_sorted(v->initial_points);
v->initial_points.clear();
if(v->l) init(init, v->l);
if(v->r) init(init, v->r);
};
init(init, root);
}
void update(Idx x, Idx y, Val z){
update_inner(x, y, z);
}
Val query(Idx lx, Idx rx, Idx ly, Idx ry){
return query_inner(rx, ly, ry) - query_inner(lx, ly, ry);
}
};
template<typename Idx, typename Val>
struct offline_static_rectangle_add_rectangle_sum{
private:
struct bit4{
int N;
std::vector<std::array<Val, 4>> sum;
bit4(int N): N(N), sum(N + 1, {0, 0, 0, 0}){}
void update(Idx l, Idx r, int lc, int rc, Val z1, Val z2){
Val a = l * z1, b = r * z1, c = l * z2, d = r * z2;
for(int i = lc + 1; i <= N; i += (i & (-i))){
sum[i][0] -= a;
sum[i][1] += z1;
sum[i][2] -= c;
sum[i][3] += z2;
}
for(int i = rc + 1; i <= N; i += (i & (-i))){
sum[i][0] += b;
sum[i][1] -= z1;
sum[i][2] += d;
sum[i][3] -= z2;
}
}
std::pair<Val, Val> query(Idx r, int rc){
Val a = 0, b = 0, c = 0, d = 0;
for(int i = rc; i > 0; i -= (i & (-i))){
a += sum[i][0];
b += sum[i][1];
c += sum[i][2];
d += sum[i][3];
}
return {a + (b * r), c + (d * r)};
}
std::pair<Val, Val> query(Idx l, Idx r, int lc, int rc){
auto [cr, dxr] = query(r, rc);
auto [cl, dxl] = query(l, lc);
return {cr - cl, dxr - dxl};
}
};
struct Update{
Idx lx, rx, ly, ry;
Val z;
int lyc, ryc;
Update(Idx lx, Idx rx, Idx ly, Idx ry, Val z, int lyc = 0, int ryc = 0): lx(lx), rx(rx), ly(ly), ry(ry), z(z), lyc(lyc), ryc(ryc){}
};
struct Query{
Idx lx, rx, ly, ry;
int id, lyc, ryc;
Query(Idx lx, Idx rx, Idx ly, Idx ry, int id, int lyc = 0, int ryc = 0): lx(lx), rx(rx), ly(ly), ry(ry), id(id), lyc(lyc), ryc(ryc){}
};
std::vector<Query> Q;
std::vector<Update> U;
void solve(std::vector<Val> &ans){
int N = U.size(), M = Q.size();
std::vector<Idx> Y;
for(int i = 0; i < N; i++){
Y.push_back(U[i].ly);
Y.push_back(U[i].ry);
}
std::sort(Y.begin(), Y.end());
Y.erase(std::unique(Y.begin(), Y.end()), Y.end());
auto lb = [&](Idx y) -> int { return std::lower_bound(Y.begin(), Y.end(), y) - Y.begin();};
for(int i = 0; i < N; i++){
int lyc = lb(U[i].ly), ryc = lb(U[i].ry);
U[i].lyc = lyc, U[i].ryc = ryc;
U.push_back(Update(U[i].rx, 0, U[i].ly, U[i].ry, -U[i].z, lyc, ryc));
}
for(int i = 0; i < M; i++){
int lyc = lb(Q[i].ly), ryc = lb(Q[i].ry);
Q[i].lyc = lyc, Q[i].ryc = ryc;
Q.push_back(Query(Q[i].rx, 0, Q[i].ly, Q[i].ry, Q[i].id + M, lyc, ryc));
}
std::sort(U.begin(), U.end(), [](const Update &a, const Update &b){return a.lx < b.lx;});
std::sort(Q.begin(), Q.end(), [](const Query &a, const Query &b){return a.lx < b.lx;});
assert(U.size() == 2 * N && Q.size() == 2 * M);
bit4 bit(Y.size());
int uid = 0, qid = 0;
while(qid < 2 * M){
if(uid < 2 * N && U[uid].lx < Q[qid].lx){
bit.update(U[uid].ly, U[uid].ry, U[uid].lyc, U[uid].ryc, -U[uid].z * U[uid].lx, U[uid].z);
uid++;
}else{
auto [a, b] = bit.query(Q[qid].ly, Q[qid].ry, Q[qid].lyc, Q[qid].ryc);
int id = Q[qid].id;
if(id >= M){
ans[id - M] += a + Q[qid].lx * b;
}else{
ans[id ] -= a + Q[qid].lx * b;
}
qid++;
}
}
}
public:
offline_static_rectangle_add_rectangle_sum(){}
void update(Idx lx, Idx rx, Idx ly, Idx ry, Val z){
U.push_back(Update(lx, rx, ly, ry, z));
}
void query(Idx lx, Idx rx, Idx ly, Idx ry){
Q.push_back(Query(lx, rx, ly, ry, Q.size()));
}
std::vector<Val> solve(){
std::vector<Val> ans(Q.size(), 0);
solve(ans);
return ans;
}
};
template<typename Val>
struct offline_dynamic_rectangle_add_rectangle_sum{
private:
struct query_type{
int lx, rx, ly, ry;
Val z;
int type, lyc, ryc;
query_type(int _lx, int _rx, int _ly, int _ry, Val _z, int _type, int _lyc = 0, int _ryc = 0):
lx(_lx), rx(_rx), ly(_ly), ry(_ry), z(_z), type(_type), lyc(_lyc), ryc(_ryc){}
};
int q = 0;
std::vector<int> qcnt{0};
std::vector<query_type> Q;
void solve(int l, int r, std::vector<Val> &ans){
if(r - l < 2) return;
int mid = (l + r) >> 1;
int left_add = (mid - l) - (qcnt[mid] - qcnt[l]);
int right_query = qcnt[r] - qcnt[mid];
if(left_add) solve(l, mid, ans);
if(right_query) solve(mid, r, ans);
if(!left_add || !right_query) return;
std::vector<query_type> Q_tmp;
// do naive
if(left_add <= 6 || right_query <= 6){
for(int i = l; i < mid; i++) if(Q[i].type == -1) Q_tmp.push_back(Q[i]);
for(int i = mid; i < r; i++){
if(Q[i].type == -1) continue;
for(query_type qi: Q_tmp){
int lx = std::max(Q[i].lx, qi.lx);
int rx = std::min(Q[i].rx, qi.rx);
if(lx >= rx) continue;
int ly = std::max(Q[i].ly, qi.ly);
int ry = std::min(Q[i].ry, qi.ry);
if(ly >= ry) continue;
ans[Q[i].type] += qi.z * Val(rx - lx) * Val(ry - ly);
}
}
return;
}
std::vector<int> Y;
for(int i = l; i < mid; i++){
if(Q[i].type != -1) continue;
Y.push_back(Q[i].ly);
Y.push_back(Q[i].ry);
}
for(int i = mid; i < r; i++){
if(Q[i].type == -1) continue;
Y.push_back(Q[i].ly);
Y.push_back(Q[i].ry);
}
std::sort(Y.begin(), Y.end());
Y.erase(std::unique(Y.begin(), Y.end()), Y.end());
auto lb = [&](int y) -> int {
return std::lower_bound(Y.begin(), Y.end(), y) - Y.begin();
};
for(int i = l; i < mid; i++){
if(Q[i].type != -1) continue;
query_type qi = Q[i];
int lyc = lb(qi.ly), ryc = lb(qi.ry);
Q_tmp.push_back(query_type(qi.lx, 0, qi.ly, qi.ry, qi.z, -1, lyc, ryc));
Q_tmp.push_back(query_type(qi.rx, 0, qi.ly, qi.ry, qi.z, -2, lyc, ryc));
}
for(int i = mid; i < r; i++){
if(Q[i].type == -1) continue;
query_type qi = Q[i];
int lyc = lb(qi.ly), ryc = lb(qi.ry);
Q_tmp.push_back(query_type(qi.lx, 0, qi.ly, qi.ry, qi.z, qi.type, lyc, ryc));
Q_tmp.push_back(query_type(qi.rx, 0, qi.ly, qi.ry, qi.z, qi.type + q, lyc, ryc));
}
std::sort(Q_tmp.begin(), Q_tmp.end(), [](const query_type &a, const query_type &b){
if(a.lx == b.lx) return a.type > b.type;
return a.lx < b.lx;
});
binary_indexed_tree_compressed<int, Val> slope(Y.size()), intercept(Y.size());
for(query_type &qi: Q_tmp){
if(qi.type == -1){
// 傾き +z, x分切片を減らす
slope.update(qi.ly, qi.ry, qi.lyc, qi.ryc, qi.z);
intercept.update(qi.ly, qi.ry, qi.lyc, qi.ryc, Val(qi.z) * Val(-qi.lx));
}else if(qi.type == -2){
// 傾き -z, x分切片を増やす
slope.update(qi.ly, qi.ry, qi.lyc, qi.ryc, Val(-qi.z));
intercept.update(qi.ly, qi.ry, qi.lyc, qi.ryc, Val(qi.z) * Val(qi.lx));
}else if(qi.type < q){
// [0, lx) × [ly, ry)を減らす
Val a = slope.query(qi.ly, qi.ry, qi.lyc, qi.ryc);
Val b = intercept.query(qi.ly, qi.ry, qi.lyc, qi.ryc);
ans[qi.type] -= Val(a) * Val(qi.lx) + Val(b);
}else{
// [0, rx) × [ly, ry)を増やす
qi.type -= q;
Val a = slope.query(qi.ly, qi.ry, qi.lyc, qi.ryc);
Val b = intercept.query(qi.ly, qi.ry, qi.lyc, qi.ryc);
ans[qi.type] += Val(a) * Val(qi.lx) + Val(b);
}
}
}
public:
offline_dynamic_rectangle_add_rectangle_sum(){}
void update(int lx, int rx, int ly, int ry, Val z){
Q.push_back(query_type(lx, rx, ly, ry, z, -1));
qcnt.push_back(0);
}
void query(int lx, int rx, int ly, int ry){
Q.push_back(query_type(lx, rx, ly, ry, 0, q++));
qcnt.push_back(1);
}
std::vector<Val> solve(){
std::vector<Val> ans(q, 0);
for(int i = 1; i < qcnt.size(); i++) qcnt[i] += qcnt[i - 1];
solve(0, Q.size(), ans);
return ans;
}
};
template<typename Idx, typename Val>
struct offline_static_cuboid_sum{
private:
offline_point_add_rectangle_sum<Idx, Val> rect;
std::vector<std::tuple<Idx, Idx, Idx, Val>> P;
struct qst{
Idx rx, ly, ry, lz, rz;
int id;
bool add;
};
std::vector<qst> Q;
std::vector<Val> __solve(){
std::sort(P.begin(), P.end());
std::sort(Q.begin(), Q.end(), [&](qst &a, qst &b){return a.rx < b.rx;});
int j = 0;
for(auto &q : Q){
while(j < P.size() && std::get<0>(P[j]) < q.rx){
auto [x, y, z, w] = P[j++];
rect.update(y, z, w);
}
rect.query(q.ly, q.ry, q.lz, q.rz);
}
std::vector<Val> ans(Q.size() / 2);
auto s = rect.solve();
j = 0;
for(auto &q : Q){
if(q.add) ans[q.id] += s[j++];
else ans[q.id] -= s[j++];
}
return ans;
}
public:
offline_static_cuboid_sum(){}
// (x, y, z)に重みwを追加
void update(Idx x, Idx y, Idx z, Val w){
P.push_back({x, y, z, w});
}
// [lx, rx) × [ly, ry) × [lz, rz)の点の重みの和
void query(Idx lx, Idx rx, Idx ly, Idx ry, Idx lz, Idx rz){
int id = Q.size() / 2;
Q.push_back({rx, ly, ry, lz, rz, id, true});
Q.push_back({lx, ly, ry, lz, rz, id, false});
}
// O((N + Q)log^2(N + Q))
std::vector<Val> solve(){
return __solve();
}
};
template<typename Idx, typename Val>
struct offline_rectangle_add_point_get{
private:
static constexpr int qlim = 1e8;
struct Query{
Idx x, y;
int id;
};
struct Update{
Idx lx, rx, ly, ry;
Val z;
};
struct Event{
Idx x;
int lyc, ryc;
Val z;
};
std::vector<Update> U;
std::vector<Query> Q;
std::vector<std::pair<bool, int>> T;
void solve(int l, int r, std::vector<Val> &ans){
if(r - l < 2) return;
int mid = (l + r) / 2;
solve(l, mid, ans);
solve(mid, r, ans);
std::vector<Idx> Y;
for(int i = mid; i < r; i++){
if(T[i].first) continue;
int id = T[i].second;
Y.push_back(Q[id].y);
}
if(Y.empty()) return;
std::sort(Y.begin(), Y.end());
Y.erase(std::unique(Y.begin(), Y.end()), Y.end());
std::vector<Event> E;
for(int i = l; i < mid; i++){
if(!T[i].first) continue;
int id = T[i].second;
int lyc = std::lower_bound(Y.begin(), Y.end(), U[id].ly) - Y.begin();
int ryc = std::lower_bound(Y.begin(), Y.end(), U[id].ry) - Y.begin();
E.push_back(Event{U[id].lx, lyc, ryc, U[id].z});
E.push_back(Event{U[id].rx, lyc, ryc, -U[id].z});
}
for(int i = mid; i < r; i++){
if(T[i].first) continue;
int id = T[i].second;
int y = std::lower_bound(Y.begin(), Y.end(), Q[id].y) - Y.begin();
E.push_back(Event{Q[id].x, y, Q[id].id + qlim, 0});
}
std::sort(E.begin(), E.end(), [](const Event &a, const Event &b){
if(a.x == b.x) return a.ryc < b.ryc;
return a.x < b.x;
});
binary_indexed_tree<Val> bit(Y.size());
for(const Event &e : E){
if(e.ryc < qlim){
if(e.lyc < Y.size()) bit.update(e.lyc, e.z);
if(e.ryc < Y.size()) bit.update(e.ryc, -e.z);
}else{
int id = e.ryc - qlim;
ans[id] += bit.query(e.lyc + 1);
}
}
}
public:
// [lx, rx) × [ly, ry)にzを足す
void update(Idx lx, Idx rx, Idx ly, Idx ry, Val z){
T.push_back({1, U.size()});
U.push_back(Update{lx, rx, ly, ry, z});
}
// get(x, y)
void query(Idx x, Idx y){
T.push_back({0, Q.size()});
Q.push_back(Query{x, y, (int)Q.size()});
}
std::vector<Val> solve(){
std::vector<Val> ans(Q.size(), 0);
solve(0, T.size(), ans);
return ans;
}
};
int main(){
prime_sieve::init(200001);
io_init();
int n;
std::cin >> n;
vector<pair<int, int>> v;
range(i, 0, n){
int a, b;
std::cin >> a >> b;
v.push_back({a, b});
}
sort(allof(v));
vector<int> mxb(200001, 0);
range(i, 0, n){
auto [a, b] = v[i];
mxb[a] = max(mxb[a], b);
}
int lower = 0;
range(i, 1, n){
if(v[i].first == v[i - 1].first){
lower = max(lower, min(v[i].second, v[i - 1].second));
}
}
// x | ai
// x | aj
// ci := ai / g
// cj := aj / g
// min(floor(bi / cj), floor(bj / ci))
// gを小さく取ってしまう -> cが大きくなる -> 答えが小さくなるため問題ない
// S[i] := iの倍数の集合
vector<vector<pair<int, int>>> S(200001);
range(i, 0, n){
auto [a, b] = v[i];
for(int d : prime_sieve::divisor(v[i].first)){
S[d].push_back({v[i].first / d, v[i].second});
if(a != d) lower = max(lower, min(b, mxb[d] / (a / d)));
}
}
// g = iと仮定していい
// {ci, bi}
// ciの条件: bjがx*ci以上
// biの条件: cjがbi/x以下
int l = lower, r = 1000000001;
while(r - l > 1){
int x = (l + r) / 2;
bool f = false;
range(i, 1, 200001){
if(f) break;
if(S[i].size() <= 10){
int sz = S[i].size();
range(j, 0, sz){
auto [ci, bi] = S[i][j];
range(k, j + 1, sz){
auto [cj, bj] = S[i][k];
if(bi >= (ll)x * cj && bj >= (ll)x * ci){
f = true;
break;
}
}
}
}else{
offline_static_rectangle_sum<ll, int> rect;
for(auto [c, b] : S[i]){
rect.update(c, b, 1);
rect.query(0, b / x + 1, (ll)x * c, 1LL << 50);
}
auto ans = rect.solve();
range(j, 0, S[i].size()){
int y = ans[j];
auto [c, b] = S[i][j];
y -= (b >= (ll)c * x ? 1 : 0);
if(y){
f = true;
break;
}
}
}
}
if(f) l = x;
else r = x;
}
std::cout << l << '\n';
}
tonegawa