D - 「火鼠的皮衣 -不焦躁的内心-」
给定
n
,
a
,
b
,
p
n,a,b,p
n,a,b,p,求下面的式子;
∑
i
=
0
⌊
n
2
⌋
a
i
b
n
−
2
i
(
n
2
i
)
m
o
d
p
\sum _{i=0}^{\lfloor {n\over 2}\rfloor } a^i b^{n-2i}{n\choose 2i} \mod p
i=0∑⌊2n⌋aibn−2i(2in)modp
数据组数 T ≤ 10000 T\le 10000 T≤10000, 1 ≤ n , a , b , p ≤ 1 0 18 1\le n,a,b,p\le 10^{18} 1≤n,a,b,p≤1018,不保证 p p p是质数。
Sol
设
a
′
=
a
(
m
o
d
p
)
a' = \sqrt a \pmod p
a′=a(modp),则原式可以化成:
∑
i
=
0
n
a
′
i
b
n
−
i
(
n
i
)
[
i is even
]
\sum_{i=0}^{n} a'^i b^{n-i} {n\choose i}[\text{i is even}]
i=0∑na′ibn−i(in)[i is even]
也就是 ( a ′ + b ) n (a'+b)^n (a′+b)n的展开式中, a ′ a' a′的次数为偶数的项的和。
设 f i , 0 / 1 f_{i,0/1} fi,0/1表示 ( a ′ + b ) i (a'+b)^i (a′+b)i中 a ′ a' a′次数为偶数/奇数的项的和,那么 f i + 1 , d = f i , d ⋅ b + f i , 1 − d ⋅ a ′ f_{i+1,d}= f_{i,d} \cdot b + f_{i,1-d} \cdot a' fi+1,d=fi,d⋅b+fi,1−d⋅a′。
矩阵快速幂就可以求出 f n , 0 f_{n,0} fn,0。
由于 a a a在模 p p p的意义下可能不存在二次剩余。但是我们在运算中会用到的数一定都可以表示成 x + a ′ ⋅ y x + a' \cdot y x+a′⋅y的形式,记录下 ( x , y ) (x,y) (x,y)就可以进行加、减、乘运算。最后得到的答案 y y y一定等于 0 0 0,输出 x x x就可以了。
还有一个问题是,这里如果用long double的那种 O ( 1 ) O(1) O(1)快速乘精度会炸, O ( log p ) O(\log p) O(logp)的慢速乘会TLE,要用int128来做乘法。
#include <cstdio>
#include <iostream>
#include <algorithm>
#include <cstring>
#define ll long long
using namespace std;
template <class T>
inline void rd(T &x) {
x=0; char c=getchar(); int f=1;
while(!isdigit(c)) { if(c=='-') f=-1; c=getchar(); }
while(isdigit(c)) x=x*10-'0'+c,c=getchar(); x*=f;
}
ll mod,W,n,a,b;
ll add(ll x) { return x>=mod?x-mod:x; }
ll dec(ll x) { return x<0?x+mod:x; }
ll Mul(ll x,ll y) {
return (__int128_t)x*y%mod;
// if(!y||(1ll<<62)/y>=x) return x*y%mod;
// ll tmp=(long double)x*y/mod;
// return ((x*y-tmp*mod)%mod+mod)%mod;
// ll res=0;
// while(y) {
// if(y&1) res=add(res+x);
// x=add(x+x),y>>=1;
// }
// return res;
}
struct Num {
ll x,y;
Num(ll x=0,ll y=0): x(x),y(y) {}
friend Num operator +(Num A,Num B) { return Num(add(A.x+B.x),add(A.y+B.y)); }
friend Num operator -(Num A,Num B) { return Num(dec(A.x-B.x),dec(A.y-B.y)); }
friend Num operator *(Num A,Num B) { return Num(add(Mul(A.x,B.x)+Mul(W,Mul(B.y,A.y))),add(Mul(A.x,B.y)+Mul(A.y,B.x))); }
};
struct Mat {
Num a[2][2];
Mat() { memset(a,0,sizeof(a)); }
void print() { for(int i=0;i<2;++i,puts("")) for(int j=0;j<2;++j) printf("(%lld,%lld) ",a[i][j].x,a[i][j].y); }
};
Mat MulMat(Mat A,Mat B) {
Mat C;
for(int i=0;i<2;++i)
for(int k=0;k<2;++k) if(A.a[i][k].x||A.a[i][k].y)
for(int j=0;j<2;++j)
C.a[i][j]=C.a[i][j]+(A.a[i][k]*B.a[k][j]);
return C;
}
Mat Pow(Mat x,ll y) {
Mat res; for(int i=0;i<2;++i) res.a[i][i]=Num(1,0);
while(y) {
if(y&1) res=MulMat(res,x);
x=MulMat(x,x),y>>=1;
}
return res;
}
int main() {
int T; rd(T);
while(T--) {
rd(n),rd(a),rd(b),rd(mod);
W=a;
Mat M;
for(int i=0;i<2;++i) for(int j=0;j<2;++j) M.a[i][j]=i==j?Num(b,0):Num(0,1);
M=Pow(M,n);
printf("%lld\n",M.a[0][0].x);
}
return 0;
}
E - 「燕的子安贝 -永命线-」
平面上有 n n n个点。你现在要画一条直线,使得到这条直线的距离不超过 d d d的点尽可能多。
输出可能的最多的点数。
n ≤ 2000 n\le 2000 n≤2000,坐标的绝对值不超过 10000 10000 10000, 0 ≤ d ≤ 10000 0\le d\le 10000 0≤d≤10000
Sol
转化一下问题,以每个点为圆心,画半径为 d d d的圆,问画一条直线至多能够交多少个圆。
那么画的直线一定与某个圆相切。枚举直线与哪个圆相切,考虑直线绕这个圆周旋转且始终与圆相切的过程中,其它的圆与这个直线相交的时间是连续的区间,求出被区间覆盖次数最大的点就可以了。
我的代码写丑了,精度有点炸,特判了 d = 0 d=0 d=0才过的。
时间复杂度 O ( n 2 log n ) O(n^2 \log n) O(n2logn)。
#include <cstdio>
#include <iostream>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <vector>
#define PB push_back
#define ll long long
#define db long double
using namespace std;
template <class T>
inline void rd(T &x) {
x=0; char c=getchar(); int f=1;
while(!isdigit(c)) { if(c=='-') f=-1; c=getchar(); }
while(isdigit(c)) x=x*10-'0'+c,c=getchar(); x*=f;
}
const int N=2010;
db eps=1e-12,Pi=acos(-1.0);
int sgn(db x) { return x<-eps?-1:(x>eps); }
struct Point {
db x,y;
Point(db x=0,db y=0): x(x),y(y) {}
friend Point operator +(Point A,Point B) { return Point(A.x+B.x,A.y+B.y); }
friend Point operator -(Point A,Point B) { return Point(A.x-B.x,A.y-B.y); }
friend Point operator *(Point A,db B) { return Point(A.x*B,A.y*B); }
}p[N];
int n,d;
struct item {
db ag; int t;
item(db ag=0,int t=0): ag(ag),t(t) {};
friend bool operator <(item A,item B) { return sgn(A.ag-B.ag)?A.ag<B.ag:A.t>B.t; }
}que[N*8];
int tot;
void add(db l,db r,int t) {
if(sgn(r-l)>=0) que[++tot]=item(l,t),que[++tot]=item(r,-t);
else add(l,Pi,t),add(-Pi,r,t);
}
db dis2(Point A) { return A.x*A.x+A.y*A.y; }
vector<db> vec;
void sol() {
int ans=0;
for(int i=1;i<=n;++i) {
vec.clear();
int cnt0=0;
for(int j=1;j<=n;++j) if(i!=j) {
if(sgn(p[i].x-p[j].x)==0) cnt0++;
else vec.PB((p[j].y-p[i].y)/(p[j].x-p[i].x));
}
ans=max(ans,cnt0+1);
sort(vec.begin(),vec.end());
for(int l=0,r;l<vec.size();l=r+1) {
r=l;
while(r+1<vec.size()&&sgn(vec[r+1]-vec[l])==0) r++;
ans=max(ans,r-l+2);
}
}
printf("%d",ans);
}
int main() {
int ans=0;
rd(n),rd(d);
for(int i=1;i<=n;++i) cin>>p[i].x>>p[i].y;
if(d==0) {
sol();
return 0;
}
for(int i=1;i<=n;++i) {
tot=0;
for(int j=1;j<=n;++j) if(j!=i) {
Point A=p[i],B=p[j];
db l=atan2(-(B.x-A.x),B.y-A.y),r=l+Pi; if(r>Pi) r-=2*Pi;
if(sgn(dis2(A-B)-4ll*d*d)>0) {
db d1=sqrt(dis2((B-A)*0.5));
db ag0=atan2(B.y-A.y,B.x-A.x);
db delt=acos(d/d1);
db a1=ag0-delt; if(-Pi>a1) a1+=2*Pi;
db a2=ag0+delt; if(a2>Pi) a2-=2*Pi;
add(l,a1,1),add(a2,r,1);
}
else add(l,r,1);
}
sort(que+1,que+tot+1);
int cnt=0;
for(int j=1;j<=tot;++j) cnt+=que[j].t,ans=max(ans,cnt+1);
}
printf("%d",ans);
return 0;
}
F - 「蓬莱的弹枝 -七色的弹幕-」
你要维护一个序列 a a a,支持以下操作:
1)将区间
[
l
,
r
]
[l,r]
[l,r]内的数向左偏移,具体来说就是
a
i
=
a
i
+
1
(
l
≤
i
<
r
)
a_i = a_{i+1} (l\le i< r )
ai=ai+1(l≤i<r),
a
r
=
a
l
a_r = a_l
ar=al。
2)将区间
[
l
,
r
]
[l,r]
[l,r]内的数全部
+
1
+1
+1。
3)给出一个下标
x
x
x,查询与
a
x
a_x
ax值相同的元素的下标与
x
x
x的差的绝对值的最小值。
n ≤ 1 0 5 , q ≤ 1 0 5 n\le 10^5,q\le 10^5 n≤105,q≤105,初始 a i ≤ 1 0 5 a_i \le 10^5 ai≤105。
Sol
对序列分块,对每个块维护:
- 加标记 t a g tag tag。
- 每个元素的真实值-加标记的值 a [ ] a[\ ] a[ ]。
- 真实值-加标记的每个值的出现次数 f [ ] f[\ ] f[ ]。
- 下标偏移值:对这个块记一个下标 s s s,表示序列中这一块的元素从左到右依次是 a [ s ] , a [ s + 1 ] ⋯ a [ m ] , a [ 1 ] , ⋯ a [ s − 1 ] a[s],a[s+1]\cdots a[m],a[1],\cdots a[s-1] a[s],a[s+1]⋯a[m],a[1],⋯a[s−1],其中 m m m是块大小。
区间加的时候这些东西都很好维护。
区间偏移的时候,对 l , r l,r l,r所属的块进行暴力,而对于中间的每一个块,相当于是删掉了最左边的元素,将其元素的下标都 − 1 -1 −1,然后在最右边插入一个元素,可以通过维护 s s s在 O ( 1 ) O(1) O(1)的时间内更新。
查询的时候,先找出左边/右边第一个出现了 a x a_x ax的块,然后在块内暴力。
时间复杂度 O ( n n ) O(n\sqrt n) O(nn)。
#include <cstdio>
#include <iostream>
#include <algorithm>
#include <cstring>
#define ll long long
using namespace std;
template <class T>
inline void rd(T &x) {
x=0; char c=getchar(); int f=1;
while(!isdigit(c)) { if(c=='-') f=-1; c=getchar(); }
while(isdigit(c)) x=x*10-'0'+c,c=getchar(); x*=f;
}
const int N=100010,M=500,C=100000;
int n;
struct BLOCK {
int a[M+10],f[N*3];
int tg,s;
int m;
void init(int *b,int _m) {
m=_m;
for(int i=0;i<m;++i) f[C+(a[i]=b[i])]++;
tg=0,s=0;
}
void add(int l,int r) { // [l, r)
if(l==0&&r==m) return (void)(tg++);
for(int i=(s+l)%m;i!=(s+r)%m;i=(i+1)%m)
f[C+a[i]]--,f[C+(++a[i])]++;
}
int flip(int y) {
int x=a[s]+tg;
f[C+a[s]]--;
f[C+(a[s]=y-tg)]++;
s=(s+1)%m;
return x;
}
int flipl(int l,int r,int x) {
int c=(s+l)%m,R=(s+r)%m,tmp=a[c]+tg;
f[C+a[c]]--;
while(c!=R) a[c]=a[(c+1)%m],c=(c+1)%m;
f[C+(a[c]=x-tg)]++;
return tmp;
}
void flipi(int l,int r) {
int tmp=a[(s+l)%m],i,R=(s+r)%m;
for(i=(s+l)%m;i!=R;i=(i+1)%m) a[i]=a[(i+1)%m];
a[i]=tmp;
}
int Qval(int p) { return a[(s+p)%m]+tg; }
int Qex(int x) { return f[C+(x-tg)]; }
int Qex_posr(int x) {
for(int i=s;;i=(i+1)%m) if(a[i]+tg==x) return (i-s+m)%m;
}
int Qex_posl(int x) {
for(int i=(s-1+m)%m;;i=(i-1+m)%m) if(a[i]+tg==x) return (i-s+m)%m;
}
int Qex_pos(int x,int l,int r,int d) {
if(l==r) return -1;
for(int i=(s+l+m)%m;i!=(s+r+m)%m;i=(i+d+m)%m) if(a[i]+tg==x) return (i-s+m)%m;
return -1;
}
}A[210];
int a[N],b[N],s[N];
void flip(int l,int r) {
if(b[l]==b[r]) return A[b[l]].flipi(l-s[b[l]],r-s[b[r]]);
int tmp=A[b[l]].Qval(l-s[b[l]]);
tmp=A[b[r]].flipl(0,r-s[b[r]],tmp);
for(int i=b[r]-1;i>b[l];--i) tmp=A[i].flip(tmp);
A[b[l]].flipl(l-s[b[l]],M-1,tmp);
}
void add(int l,int r) {
if(b[l]==b[r]) return A[b[l]].add(l-s[b[l]],r-s[b[r]]+1);
A[b[l]].add(l-s[b[l]],M);
A[b[r]].add(0,r-s[b[r]]+1);
for(int i=b[l]+1;i<b[r];++i) A[i].tg++;
}
int query(int x) {
int v=A[b[x]].Qval(x-s[b[x]]);
int R=A[b[x]].Qex_pos(v,x-s[b[x]]+1,0,1);
int L=A[b[x]].Qex_pos(v,x-s[b[x]]-1,-1,-1);
if(L==-1) {
int p=b[x]-1;
while(p>=0&&!A[p].Qex(v)) p--;
if(p>=0) L=s[p]+A[p].Qex_posl(v);
else L=-1e8;
}
else L+=s[b[x]];
if(R==-1) {
int p=b[x]+1;
while(p<=b[n-1]&&!A[p].Qex(v)) p++;
if(p<=b[n-1]) R=s[p]+A[p].Qex_posr(v);
else R=1e8;
}
else R+=s[b[x]];
int ans=min(x-L,R-x);
return ans<=n?ans:-1;
}
int main() {
int q; rd(n),rd(q);
for(int i=0;i<n;++i) rd(a[i]),b[i]=i/M;
for(int p=0;p<n;p+=M) A[p/M].init(a+p,min(M,n-p)),s[p/M]=(p/M)*M;
while(q--) {
int op,l,r; rd(op);
if(op==1) rd(l),rd(r),flip(--l,--r);
else if(op==2) rd(l),rd(r),add(--l,--r);
else rd(l),printf("%d\n",query(--l));
}
return 0;
}
ans<=n?ans:-1;
}
int main() {
int q; rd(n),rd(q);
for(int i=0;i<n;++i) rd(a[i]),b[i]=i/M;
for(int p=0;p<n;p+=M) A[p/M].init(a+p,min(M,n-p)),s[p/M]=(p/M)*M;
while(q--) {
int op,l,r; rd(op);
if(op==1) rd(l),rd(r),flip(--l,--r);
else if(op==2) rd(l),rd(r),add(--l,--r);
else rd(l),printf("%d\n",query(--l));
}
return 0;
}