这题思路很清楚。
应该就是先把b数组搞出来,然后再根据b算c数组的中位数。
首先先把b搞出来,如果用倍增,那么预处理(O(n log n log C)),算出b复杂度为(O(n log n log C log C)),稍微用一下gcd的trick就可以变成预处理(O(n (log n+log C))) ,算b复杂度(O(n log n log C))
如果用线段树加trick维护 预处理(O(n log C)),算b在树上二分也是(O(n log n log C))的。
然而我用了裸的线段树,所以是三个log的。
然后比较恶心的是类欧的部分。
假设要求
(sum_{i=0}^n lfloor frac{ai+b}{c}
floor)
(sum_{i=0}^n sum_{j=1} [lfloor frac{ai+b}{c}
floor geq j])
(sum_{i=0}^n sum_{j=1} [frac{ai+b}{c} geq j])
(sum_{i=0}^n sum_{j=0} [frac{ai+b}{c} geq j+1])
(sum_{i=0}^n sum_{j=0} [ai geq j*c+c-b])
(sum_{j=0} sum_{i=0}^n [i geq frac{j*c+c-b}{a}])
(sum_{j=0} sum_{i=0}^n [i > frac{j*c+c-b-1}{a}])
(sum_{j=0} n- lfloor frac{j*c+c-b-1}{a}
floor)
令(m=lfloor frac{an+b}{c}
floor)
(n*m-sum_{j=0} lfloor frac{j*c+c-b-1}{a}
floor)
#include <bits/stdc++.h>
using namespace std;
typedef pair<int,int> p;
typedef long long ll;
namespace Segment{
//also bz
const int P=140000;
int T[P];
void update(int x,int l,int r,int pos,int v){
if (l==r){
T[x]=v;
return;
}
int mid=(l+r)>>1;
if (pos<=mid) update(x<<1,l,mid,pos,v);
else update(x<<1|1,mid+1,r,pos,v);
T[x]=__gcd(T[x<<1],T[x<<1|1]);
}
int ask(int x,int l,int r,int ll,int rr){
if (ll<=l&&r<=rr) return T[x];
int mid=(l+r)>>1,ret=0;
if (ll<=mid) ret=__gcd(ret,ask(x<<1,l,mid,ll,rr));
if (mid<rr) ret=__gcd(ret,ask(x<<1|1,mid+1,r,ll,rr));
return ret;
}
}
namespace Math{
ll stand(ll lim1,ll a,ll b,ll c){
//cerr<<"stand"<<lim1<<" "<<a<<" "<<b<<" "<<c<<endl;
// x[0,lim1] ax+b/c
if (a<0){
ll pl=-(((a+1)/c)-1);
return stand(lim1,a+pl*c,b,c)-pl*lim1*(lim1+1)/2;
}
if (a==0) return (lim1+1)*(b/c);
if (a>=c||b>=c){
return stand(lim1,a%c,b%c,c)+(a/c)*lim1*(lim1+1)/2+(b/c)*(lim1+1);
}
ll m=(a*lim1+b)/c;
return lim1*m-stand(m-1,c,c-b-1,a);
}
ll divide(ll x,ll y){
if (x>=0) return x/y;
return x%y?x/y-1:x/y;
}
ll oc(ll a,ll b,ll lim1,ll lim2,ll c){
//cerr<<"oc"<<a<<" "<<b<<" "<<lim1<<" "<<lim2<<" "<<c<<endl;
if (lim1<0||lim2<0||c<0) return 0;
//ax+by<=c x [0,lim1] y [0,lim2]
lim1=min(lim1,c/a);//max 0
ll ret=lim1+1;//y=0
//cerr<<"ret"<<ret<<endl;
ll cd=min(max(divide(c-b*lim2,a),-1ll),lim1);
/*for (int l=0,r=lim1,mid=(l+r)>>1; l<=r; mid=(l+r)>>1)
if ((c-a*mid)/b>=lim2) cd=mid,l=mid+1; else r=mid-1;*/
ret+=(cd+1)*lim2;
//cout<<"ret"<<ret<<" "<<cd<<endl;
ret+=stand(lim1,-a,c,b);
if (cd>-1) ret-=stand(cd,-a,c,b);
//cerr<<"ret"<<ret<<endl;
return ret;
}
ll co(ll a,ll b,ll c){
//cerr<<"CP"<<a<<" "<<b<<" "<<c<<endl;
ll num=b/a;
if (num>c) return c*(c+1)/2;
return (c-num)*num+num*(num+1)/2;
}
}
const int N=50010;
const int C=100000;
ll b[C+10],sum[C+10],sb[C+10];
int pre[C+10];
ll calc(ll lim){
//cerr<<"calc"<<lim<<endl;
ll ret=0;
for (int i=1; i<=C; ++i) if (b[i]) ret+=Math::co(i,lim-1,b[i]);
for (int l=1,r=1; l<C; ++l){
if (!b[l]) continue;
if (l<r) ret+=b[l]*(sb[r-1]-sb[l]);
for (; r<=C&&sum[r]-sum[l]<lim; ++r)
if (l!=r&&b[r]) ret+=Math::oc(l,r,b[l]-1,b[r]-1,lim-1-(sum[r-1]-sum[l])-l-r);
if (r<=C){
if (l!=r&&b[r]) ret+=Math::oc(l,r,b[l]-1,b[r]-1,lim-1-(sum[r-1]-sum[l])-l-r);
}
}
return ret;
}
int n;
int main(){
//cerr<<Math::oc(2,5,0,1,4)<<endl;
//return 0;
scanf("%d",&n);
for (int i=1,x; i<=n; ++i){
scanf("%d",&x);
int la=i;
Segment::update(1,1,n,i,x);
while (la){
int tmp=la;
for (int l=1,r=la,mid=(l+r)>>1; l<=r; mid=(l+r)>>1)
if (Segment::ask(1,1,n,mid,i)==Segment::ask(1,1,n,la,i)) la=mid,r=mid-1; else l=mid+1;
b[Segment::ask(1,1,n,la,i)]+=tmp-la+1;
--la;
}
}
//for (int i=1; i<=20; ++i) cout<<b[i]<<" "; cout<<endl;
for (int i=1; i<=C; ++i){
sb[i]=sb[i-1]+b[i];
sum[i]=sum[i-1]+b[i]*i;
}
int la=0;
for (int i=1; i<=C; ++i)
if (b[i]){
pre[i]=la;
la=i;
}
ll bat=(ll)n*(n+1)/2,ans;
bat=bat*(bat+1)/2;
//cerr<<"BAT"<<bat<<endl;
//cerr<<calc(8)<<endl;
//return 0;
for (ll l=1,r=(ll)C*n*(n+1)/2,mid=(l+r)>>1; l<=r; mid=(l+r)>>1){
//cerr<<"NA"<<ans<<endl;
if (calc(mid)+1<=(bat+1)/2) ans=mid,l=mid+1;
else r=mid-1;
}
cout<<ans;
}