#include <bits/stdc++.h>
using namespace std;
struct sum_kth_smallest {
struct Node {
long long sum;
int cnt;
int lCh, rCh; //children, indexes into `tree`
} ;
int mn, mx;
vector< int > roots;
deque< Node> tree;
sum_kth_smallest( const vector< int > & arr) : mn( INT_MAX ) , mx( INT_MIN ) , roots( arr.size ( ) + 1 , 0 ) {
tree.push_back ( { 0 , 0 , 0 } ) ; //acts as null
for ( int val : arr) mn = min( mn, val) , mx = max( mx, val) ;
for ( int i = 0 ; i < ( int ) arr.size ( ) ; i++ )
roots[ i + 1 ] = update( roots[ i] , - mx, mx, arr[ i] ) ;
}
int update( int v, int tl, int tr, int idx) {
if ( tl == tr) {
tree.push_back ( { tree[ v] .sum + tl, tree[ v] .cnt + 1 , 0 , 0 } ) ;
return tree.size ( ) - 1 ;
}
int tm = tl + ( tr - tl) / 2 ;
int lCh = tree[ v] .lCh ;
int rCh = tree[ v] .rCh ;
if ( idx <= tm )
lCh = update( lCh, tl, tm , idx) ;
else
rCh = update( rCh, tm + 1 , tr, idx) ;
tree.push_back ( { tree[ lCh] .sum + tree[ rCh] .sum , tree[ lCh] .cnt + tree[ rCh] .cnt , lCh, rCh} ) ;
return tree.size ( ) - 1 ;
}
/* find kth smallest number among arr[l], arr[l+1], ..., arr[r]
* k is 1-based, so find_kth(l,r,1) returns the min
*/
int query( int l, int r, int k) const {
assert ( 1 <= k && k <= r - l + 1 ) ; //note this condition implies L <= R
assert ( 0 <= l && r + 1 < ( int ) roots.size ( ) ) ;
return query( roots[ l] , roots[ r + 1 ] , - mx, mx, k) ;
}
int query( int vl, int vr, int tl, int tr, int k) const {
if ( tl == tr)
return tl;
int tm = tl + ( tr - tl) / 2 ;
int left_count = tree[ tree[ vr] .lCh ] .cnt - tree[ tree[ vl] .lCh ] .cnt ;
if ( left_count >= k) return query( tree[ vl] .lCh , tree[ vr] .lCh , tl, tm , k) ;
return query( tree[ vl] .rCh , tree[ vr] .rCh , tm + 1 , tr, k - left_count) ;
}
/* find **sum** of k smallest numbers among arr[l], arr[l+1], ..., arr[r]
* k is 1-based, so find_kth(l,r,1) returns the min
*/
long long query_sum( int l, int r, int k) const {
assert ( 1 <= k && k <= r - l + 1 ) ; //note this condition implies L <= R
assert ( 0 <= l && r + 1 < ( int ) roots.size ( ) ) ;
return query_sum( roots[ l] , roots[ r + 1 ] , - mx, mx, k) ;
}
long long query_sum( int vl, int vr, int tl, int tr, int k) const {
if ( tl == tr)
return 1LL * tl * k;
int tm = tl + ( tr - tl) / 2 ;
int left_count = tree[ tree[ vr] .lCh ] .cnt - tree[ tree[ vl] .lCh ] .cnt ;
long long left_sum = tree[ tree[ vr] .lCh ] .sum - tree[ tree[ vl] .lCh ] .sum ;
if ( left_count >= k) return query_sum( tree[ vl] .lCh , tree[ vr] .lCh , tl, tm , k) ;
return left_sum + query_sum( tree[ vl] .rCh , tree[ vr] .rCh , tm + 1 , tr, k - left_count) ;
}
} ;
//MUCH RANDOM!!!
seed_seq seed{
( uint32_t ) chrono:: duration_cast < chrono:: nanoseconds > ( chrono:: high_resolution_clock :: now ( ) .time_since_epoch ( ) ) .count ( ) ,
( uint32_t ) random_device( ) ( ) ,
( uint32_t ) ( uintptr_t ) make_unique< char > ( ) .get ( ) ,
( uint32_t ) __builtin_ia32_rdtsc( )
} ;
mt19937 rng( seed) ;
template < class T>
inline T getRand( T l, T r) {
assert ( l <= r) ;
return uniform_int_distribution< T> ( l, r) ( rng) ;
}
int main( ) {
while ( true ) {
int n = getRand( 1 , 1000 ) ;
cout << "start of new test. n = " << n << endl;
vector< int > arr( n) ;
for ( int i = 0 ; i < n; i++ ) {
arr[ i] = getRand< int > ( 0 , 1e5 ) ;
}
sum_kth_smallest st( arr) ;
for ( int queries = 1000 ; queries-- ; ) {
int L = getRand( 0 ,n- 1 ) , R = getRand( 0 ,n- 1 ) ;
if ( L > R) swap( L,R) ;
vector< int > subarr( R- L+ 1 ) ;
copy( arr.begin ( ) + L, arr.begin ( ) + R+ 1 , subarr.begin ( ) ) ;
sort( subarr.begin ( ) , subarr.end ( ) ) ;
int numLess = 0 ;
long long prefixSum = 0 ;
for ( int k = 1 ; k <= R- L+ 1 ; k++ ) {
prefixSum + = subarr[ k- 1 ] ;
assert ( st.query ( L,R,k) == subarr[ k- 1 ] ) ;
assert ( st.query_sum ( L,R,k) == prefixSum) ;
}
}
}
return 0 ;
}
I2luY2x1ZGUgPGJpdHMvc3RkYysrLmg+CnVzaW5nIG5hbWVzcGFjZSBzdGQ7CgoKc3RydWN0IHN1bV9rdGhfc21hbGxlc3QgewoKCXN0cnVjdCBOb2RlIHsKCQlsb25nIGxvbmcgc3VtOwoJCWludCBjbnQ7CgkJaW50IGxDaCwgckNoOy8vY2hpbGRyZW4sIGluZGV4ZXMgaW50byBgdHJlZWAKCX07CgoJaW50IG1uLCBteDsKCXZlY3RvcjxpbnQ+IHJvb3RzOwoJZGVxdWU8Tm9kZT4gdHJlZTsKCglzdW1fa3RoX3NtYWxsZXN0KGNvbnN0IHZlY3RvcjxpbnQ+JiBhcnIpIDogbW4oSU5UX01BWCksIG14KElOVF9NSU4pLCByb290cyhhcnIuc2l6ZSgpICsgMSwgMCkgewoJCXRyZWUucHVzaF9iYWNrKHswLCAwLCAwfSk7IC8vYWN0cyBhcyBudWxsCgkJZm9yIChpbnQgdmFsIDogYXJyKSBtbiA9IG1pbihtbiwgdmFsKSwgbXggPSBtYXgobXgsIHZhbCk7CgkJZm9yIChpbnQgaSA9IDA7IGkgPCAoaW50KWFyci5zaXplKCk7IGkrKykKCQkJcm9vdHNbaSArIDFdID0gdXBkYXRlKHJvb3RzW2ldLCAtbXgsIG14LCBhcnJbaV0pOwoJfQoJaW50IHVwZGF0ZShpbnQgdiwgaW50IHRsLCBpbnQgdHIsIGludCBpZHgpIHsKCQlpZiAodGwgPT0gdHIpIHsKCQkJdHJlZS5wdXNoX2JhY2soe3RyZWVbdl0uc3VtICsgdGwsIHRyZWVbdl0uY250ICsgMSwgMCwgMH0pOwoJCQlyZXR1cm4gdHJlZS5zaXplKCkgLSAxOwoJCX0KCQlpbnQgdG0gPSB0bCArICh0ciAtIHRsKSAvIDI7CgkJaW50IGxDaCA9IHRyZWVbdl0ubENoOwoJCWludCByQ2ggPSB0cmVlW3ZdLnJDaDsKCQlpZiAoaWR4IDw9IHRtKQoJCQlsQ2ggPSB1cGRhdGUobENoLCB0bCwgdG0sIGlkeCk7CgkJZWxzZQoJCQlyQ2ggPSB1cGRhdGUockNoLCB0bSArIDEsIHRyLCBpZHgpOwoJCXRyZWUucHVzaF9iYWNrKHt0cmVlW2xDaF0uc3VtICsgdHJlZVtyQ2hdLnN1bSwgdHJlZVtsQ2hdLmNudCArIHRyZWVbckNoXS5jbnQsIGxDaCwgckNofSk7CgkJcmV0dXJuIHRyZWUuc2l6ZSgpIC0gMTsKCX0KCgoJLyogZmluZCBrdGggc21hbGxlc3QgbnVtYmVyIGFtb25nIGFycltsXSwgYXJyW2wrMV0sIC4uLiwgYXJyW3JdCgkgKiBrIGlzIDEtYmFzZWQsIHNvIGZpbmRfa3RoKGwsciwxKSByZXR1cm5zIHRoZSBtaW4KCSAqLwoJaW50IHF1ZXJ5KGludCBsLCBpbnQgciwgaW50IGspIGNvbnN0IHsKCQlhc3NlcnQoMSA8PSBrICYmIGsgPD0gciAtIGwgKyAxKTsgLy9ub3RlIHRoaXMgY29uZGl0aW9uIGltcGxpZXMgTCA8PSBSCgkJYXNzZXJ0KDAgPD0gbCAmJiByICsgMSA8IChpbnQpcm9vdHMuc2l6ZSgpKTsKCQlyZXR1cm4gcXVlcnkocm9vdHNbbF0sIHJvb3RzW3IgKyAxXSwgLW14LCBteCwgayk7Cgl9CglpbnQgcXVlcnkoaW50IHZsLCBpbnQgdnIsIGludCB0bCwgaW50IHRyLCBpbnQgaykgY29uc3QgewoJCWlmICh0bCA9PSB0cikKCQkJcmV0dXJuIHRsOwoJCWludCB0bSA9IHRsICsgKHRyIC0gdGwpIC8gMjsKCQlpbnQgbGVmdF9jb3VudCA9IHRyZWVbdHJlZVt2cl0ubENoXS5jbnQgLSB0cmVlW3RyZWVbdmxdLmxDaF0uY250OwoJCWlmIChsZWZ0X2NvdW50ID49IGspIHJldHVybiBxdWVyeSh0cmVlW3ZsXS5sQ2gsIHRyZWVbdnJdLmxDaCwgdGwsIHRtLCBrKTsKCQlyZXR1cm4gcXVlcnkodHJlZVt2bF0uckNoLCB0cmVlW3ZyXS5yQ2gsIHRtICsgMSwgdHIsIGsgLSBsZWZ0X2NvdW50KTsKCX0KCgkvKiBmaW5kICoqc3VtKiogb2YgayBzbWFsbGVzdCBudW1iZXJzIGFtb25nIGFycltsXSwgYXJyW2wrMV0sIC4uLiwgYXJyW3JdCgkgKiBrIGlzIDEtYmFzZWQsIHNvIGZpbmRfa3RoKGwsciwxKSByZXR1cm5zIHRoZSBtaW4KCSAqLwoJbG9uZyBsb25nIHF1ZXJ5X3N1bShpbnQgbCwgaW50IHIsIGludCBrKSBjb25zdCB7CgkJYXNzZXJ0KDEgPD0gayAmJiBrIDw9IHIgLSBsICsgMSk7IC8vbm90ZSB0aGlzIGNvbmRpdGlvbiBpbXBsaWVzIEwgPD0gUgoJCWFzc2VydCgwIDw9IGwgJiYgciArIDEgPCAoaW50KXJvb3RzLnNpemUoKSk7CgkJcmV0dXJuIHF1ZXJ5X3N1bShyb290c1tsXSwgcm9vdHNbciArIDFdLCAtbXgsIG14LCBrKTsKCX0KCWxvbmcgbG9uZyAgcXVlcnlfc3VtKGludCB2bCwgaW50IHZyLCBpbnQgdGwsIGludCB0ciwgaW50IGspIGNvbnN0IHsKCQlpZiAodGwgPT0gdHIpCgkJCXJldHVybiAxTEwgKiB0bCAqIGs7CgkJaW50IHRtID0gdGwgKyAodHIgLSB0bCkgLyAyOwoJCWludCBsZWZ0X2NvdW50ID0gdHJlZVt0cmVlW3ZyXS5sQ2hdLmNudCAtIHRyZWVbdHJlZVt2bF0ubENoXS5jbnQ7CgkJbG9uZyBsb25nIGxlZnRfc3VtID0gdHJlZVt0cmVlW3ZyXS5sQ2hdLnN1bSAtIHRyZWVbdHJlZVt2bF0ubENoXS5zdW07CgkJaWYgKGxlZnRfY291bnQgPj0gaykgcmV0dXJuIHF1ZXJ5X3N1bSh0cmVlW3ZsXS5sQ2gsIHRyZWVbdnJdLmxDaCwgdGwsIHRtLCBrKTsKCQlyZXR1cm4gbGVmdF9zdW0gKyBxdWVyeV9zdW0odHJlZVt2bF0uckNoLCB0cmVlW3ZyXS5yQ2gsIHRtICsgMSwgdHIsIGsgLSBsZWZ0X2NvdW50KTsKCX0KfTsKCgovL01VQ0ggUkFORE9NISEhCnNlZWRfc2VxIHNlZWR7CgkodWludDMyX3QpY2hyb25vOjpkdXJhdGlvbl9jYXN0PGNocm9ubzo6bmFub3NlY29uZHM+KGNocm9ubzo6aGlnaF9yZXNvbHV0aW9uX2Nsb2NrOjpub3coKS50aW1lX3NpbmNlX2Vwb2NoKCkpLmNvdW50KCksCgkodWludDMyX3QpcmFuZG9tX2RldmljZSgpKCksCgkodWludDMyX3QpKHVpbnRwdHJfdCltYWtlX3VuaXF1ZTxjaGFyPigpLmdldCgpLAoJKHVpbnQzMl90KV9fYnVpbHRpbl9pYTMyX3JkdHNjKCkKfTsKbXQxOTkzNyBybmcoc2VlZCk7Cgp0ZW1wbGF0ZTxjbGFzcyBUPgppbmxpbmUgVCBnZXRSYW5kKFQgbCwgVCByKSB7Cglhc3NlcnQobCA8PSByKTsKCXJldHVybiB1bmlmb3JtX2ludF9kaXN0cmlidXRpb248VD4obCwgcikocm5nKTsKfQoKaW50IG1haW4oKSB7Cgl3aGlsZSh0cnVlKSB7CgkJaW50IG4gPSBnZXRSYW5kKDEsIDEwMDApOwoJCWNvdXQgPDwgInN0YXJ0IG9mIG5ldyB0ZXN0LiBuID0gIiA8PCBuIDw8IGVuZGw7CgkJdmVjdG9yPGludD4gYXJyKG4pOwoJCWZvcihpbnQgaSA9IDA7IGkgPCBuOyBpKyspIHsKCQkJYXJyW2ldID0gZ2V0UmFuZDxpbnQ+KDAsIDFlNSk7CgkJfQoJCXN1bV9rdGhfc21hbGxlc3Qgc3QoYXJyKTsKCQlmb3IoaW50IHF1ZXJpZXMgPSAxMDAwOyBxdWVyaWVzLS07KSB7CgkJCWludCBMID0gZ2V0UmFuZCgwLG4tMSksIFIgPSBnZXRSYW5kKDAsbi0xKTsKCQkJaWYoTCA+IFIpIHN3YXAoTCxSKTsKCQkJdmVjdG9yPGludD4gc3ViYXJyKFItTCsxKTsKCQkJY29weShhcnIuYmVnaW4oKStMLCBhcnIuYmVnaW4oKStSKzEsIHN1YmFyci5iZWdpbigpKTsKCQkJc29ydChzdWJhcnIuYmVnaW4oKSwgc3ViYXJyLmVuZCgpKTsKCQkJaW50IG51bUxlc3MgPSAwOwoJCQlsb25nIGxvbmcgcHJlZml4U3VtID0gMDsKCQkJZm9yKGludCBrID0gMTsgayA8PSBSLUwrMTsgaysrKSB7CgkJCQlwcmVmaXhTdW0gKz0gc3ViYXJyW2stMV07CgkJCQlhc3NlcnQoc3QucXVlcnkoTCxSLGspID09IHN1YmFycltrLTFdKTsKCQkJCWFzc2VydChzdC5xdWVyeV9zdW0oTCxSLGspID09IHByZWZpeFN1bSk7CgkJCX0KCQl9Cgl9CglyZXR1cm4gMDsKfQo=