题解:P10666 COUNTARI
题目传送门
神秘做法。
题目分析
记
首先我们先不算三个数相等的情况,最后再加上去。这样方便后面处理。
先不考虑位置顺序的限制,发现是好求的。您都开到这道题了那一定会的对吧!。记
考虑减掉多算的部分,发现只有两种情况:
记现在处理到
-
将
f 自己卷起来,得到h ,答案为h_{2a_p} 。 -
枚举前面的每个数,逐个求出答案,为
\sum f_{2a_p-a_i} 。
发现上面两种方法都是很慢的。将两种方法结合一下,类似于分块,我们可以快速将前几个块末求出的
设块长为
于是就做完啦。
Code
// by wangyizhi(571247)
#include<bits/stdc++.h>
//#include<bits/extc++.h>
using namespace std;
using ll=long long;
using ld=long double;
//#define int ll
using pii=pair<int,int>;
//bool Mst;
namespace wyzmath
{
/*complex*/
template<typename _Tp=double>
class complex
{
public:
typedef _Tp vtype;
vtype x,y;
complex(vtype _x=0,vtype _y=0):x(_x),y(_y){}
vtype length(){return std::sqrt(x*x+y*y);}
}; // class complex
/*about complex*/
template<typename _Tp=double>
complex<_Tp> operator+(complex<_Tp> a,complex<_Tp> b)
{return complex<_Tp>(a.x+b.x,a.y+b.y);}
template<typename _Tp=double>
complex<_Tp> operator-(complex<_Tp> a,complex<_Tp> b)
{return complex<_Tp>(a.x-b.x,a.y-b.y);}
template<typename _Tp=double>
complex<_Tp> operator*(complex<_Tp> a,complex<_Tp> b)
{return complex<_Tp>(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}
template<typename _Tp=double>
complex<_Tp> operator+=(complex<_Tp> &a,complex<_Tp> b){return a=a+b;}
template<typename _Tp=double>
complex<_Tp> operator-=(complex<_Tp> &a,complex<_Tp> b){return a=a-b;}
template<typename _Tp=double>
complex<_Tp> operator*=(complex<_Tp> &a,complex<_Tp> b){return a=a*b;}
}
using namespace std;
using ld=long double;
using compld=wyzmath::complex<ld>;
const compld I(1,0);
const int N=1<<21,V=3e4,VV=6e4,B=5000;
const ld PI=3.141592653589793;
int rev[N];
ll a[N],c[N],s[N],d[N],e[N],p[N];
void FFT(compld *a,int n,int op)
{
for(int i=0;i<n;i++) if(i<rev[i]) std::swap(a[i],a[rev[i]]);
for(int mid=1;mid<n;mid<<=1)
{
compld omega(cos(PI/mid),op*sin(PI/mid));
for(int len=mid<<1,i=0;i<n;i+=len)
{
compld o(1,0),x,y;
for(int j=0;j<mid;j++,o=o*omega)
x=a[i+j],y=o*a[i+j+mid],a[i+j]=x+y,a[i+j+mid]=x-y;
}
}
}
void conv(ll *a,ll *b,ll *c,int n)
{
static compld tmpa[N],tmpb[N];
int sz=1,l=0;
while(sz<n) sz<<=1,l++;
for(int i=1;i<sz;i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
for(int i=0;i<sz;i++) tmpa[i]=a[i],tmpb[i]=b[i];
FFT(tmpa,sz,1),FFT(tmpb,sz,1);
for(int i=0;i<sz;i++) tmpa[i]=tmpa[i]*tmpb[i];
FFT(tmpa,sz,-1);
for(int i=0;i<sz;i++) c[i]=(ll)(tmpa[i].x/sz+0.5);
for(int i=n;i<sz;i++) c[i]=0;
}
//bool Med;
signed main()
{
// cerr<<"Memory Size: "<<abs((&Med)-(&Mst))/1024.0/1024<<" MB\n";
ios::sync_with_stdio(0);cin.tie(0),cout.tie(0);
// freopen("in.in","r",stdin);
// freopen("out.out","w",stdout);
int n;ll ans=0;
cin>>n;
for(int i=1;i<=n;i++) cin>>a[i],a[i]--,c[a[i]]++;
conv(c,c,s,VV);
for(int i=0;i<V;i++) s[i<<1]-=c[i]*c[i];
for(int i=0;i<V;i++) ans+=s[i<<1]*c[i];
// cout<<ans<<"\n";
auto sol=[&]()->void
{
for(int i=0;i<VV;i++) s[i]=0;
for(int i=0;i<V;i++) d[i]=e[i]=0;
for(int i=1,le=0;i<=n;i++)
{
ans-=s[a[i]*2]-d[a[i]]*d[a[i]];
for(int j=le+1;j<=i;j++) if(a[i]!=a[j]&&a[i]*2>=a[j])
ans-=d[a[i]*2-a[j]]*2+e[a[i]*2-a[j]];
e[a[i]]++;
if(i%B==0&&i!=n)
{
for(int j=0;j<V;j++) d[j]+=e[j],e[j]=0;
conv(d,d,s,VV);
le=i;
}
}
};
sol();reverse(a+1,a+n+1);sol();ans/=2;
for(int i=0;i<V;i++) ans+=c[i]*(c[i]-1)*(c[i]-2)/6;
cout<<ans<<"\n";
return 0;
}