第一部分:分析题意
钻孔实际上可以看作令$b_{x}$对$h$取$min$,当$b_{x}le h$时无意义,所以不妨假设$b_{x}>h$
同时,若$max(a_{x},a_{x+1})le h$,显然也不会发生流动,所以不妨假设$max(a_{x},a_{x+1})>h$
根据对称性,$a_{x}>a_{x+1}$和$a_{x}<a_{x+1}$也是类似的,不妨仅考虑前者
简单归纳,可以发现此时流动一定是从左到右的,且发生流动的必然是一个连续区间
考虑确定流动区间$[l,r]$,使得$forall lle i<r$,$a_{i}$向$a_{i+1}$发生流动(且满足此条件下区间最大)
当确定此区间以及最终$a_{r}$的值$h$后,根据流动的性质,即可得到$forall lle ile r,a_{i}=max(a_{i+1},b_{i})$,对后者进行归纳也即为$max(h,max_{j=i}^{r-1}b_{j})$
称此操作为$ML(l,r,h)$,再记此操作所需水量为$FL(l,r,h)=sum_{i=l}^{r}max(h,max_{j=i}^{r-1}b_{j})$
下面,即来考虑确定这个区间$[l,r]$和$h$的值——
关于流动区间的左端点$l$,显然为最大的$lle x$,使得$b_{l-1}ge a_{x}$(约定$b_{0}=infty$)
(注意由于挡板比$a_{x}$小,若该位置水量不为$a_{x}$则会发生流动,因此有$a_{l}=a_{l+1}=...=a_{x}$)
关于流动区间的右端点$r$,注意到当一个区间满足$FL(l,r,a_{r})le sum_{i=l}^{r}a_{i}$是其$a_{r}$被影响的必要条件,也即要找到最大的$r$,使得满足此条件
(关于这里,更好的是取小于号,但由于是实数,这里就取小于等于了)
类似地,对于$h$也即找到最大的$h$,使得$FL(l,r,h)le sum_{i=l}^{r}a_{i}$(即恰好取到等号时)
(上面仅仅是感性的分析此问题,下面将给出比较严谨的描述)
第二部分:建立模型
综上分析,修改$b_{x}$为$h$后,即执行以下过程:
1.找到最大的$lle x$,使得$b_{l-1}ge a_{x}$(约定$b_{0}=infty$)
2.找到最大的$r$,使得$FL(l,r,a_{r})le sum_{i=l}^{r}a_{i}$
3.找到最大的$h$,使得$FL(l,r,h)le sum_{i=l}^{r}a_{i}$
4.执行操作$ML(l,r,h)$
其中$l$是最为简单的,用线段树维护$b_{i}$区间最大值即可,复杂度为$o(qlog n)$,以下不考虑(即假设已知$l$)
对于$r$和$h$,其是具有单调性的,其中$h$是显然的,下面来说明$r$——
更具体的,即求证若$FL(l,r,a_{r})le sum_{i=l}^{r}a_{i}$,则$FL(l,r-1,a_{r-1})le sum_{i=l}^{r-1}a_{i}$(以此类推,所有比$r$小的右端点都会满足)
反证法,若$FL(l,r-1,a_{r-1})>sum_{i=l}^{r-1}a_{i}$,即可推出$FL(l,r,a_{r})-FL(l,r-1,a_{r-1})<a_{r}$
同时,注意到$a_{r-1}le max(b_{r-1},a_{r})$,否则即会发生流动,因此显然矛盾
由此,对于$r$和$h$二分(后者在实数范围内),并$o(n)$求出$FL(l,r,a_{r})$和$sum_{i=l}^{r}a_{i}$即可
对于$ML(l,r,h)$,同样根据定义执行即可,复杂度为$o(nq)$
第三部分:线段树优化
考虑使用线段树来优化$FL(l,r,h)$和$ML(l,r,h)$的复杂度,具体即——
假设在线段树上递归到区间$[l,r]$(被完全包含),以$FL(l,r,h)$为例,即先递归$FL(mid+1,r,h)$,并得到新的$h'=max(h,max_{i=mid}^{r-1}b_{i})$,再递归左区间$FL(l,mid,h')$
对于线段树的每一个节点,预处理出$max_{i=mid}^{r-1}b_{i}$,每次修改再暴力在$b_{i}$的线段树上查询,显然单次至多影响$o(log n)$个,维护这个最大值的复杂度为$o(log^{2}n)$
由此,我们即可在一个位置上$o(1)$算出$h'$,并对$h'$的值分类讨论:
1.若$h'=h$,显然右区间的结果应为$(r-mid)h$
2.若$h'=max_{i=mid}^{r-1}b_{i}$,则左区间即查询$FL(l,mid,max_{i=mid}^{r-1}b_{i})$,与$h$无关
对于第2类,再在每一个节点上预处理出$FL(l,mid,max_{i=mid}^{r-1}b_{i})$的值,那么单次$FL(l,r,h)$的查询复杂度即变为$o(log^{2}n)$
(即使当前区间已经被完整覆盖,仍要递归其中一个儿子,复杂度显然即$o(log^{2}n)$)
但特别的,若$FL(l,r,h)$中$[l,r]$恰为线段树上的某个区间时,复杂度降为$o(log n)$
(另外,对于未被完全包含的区间$[l,r]$,可以先向右递归,并得到这个最大值)
下面,考虑在修改时维护$FL(l,mid,max_{i=mid}^{r-1}b_{i})$,与$max_{i=mid}^{r-1}b_{i}$相同,至多影响$o(log n)$个,且每一个计算复杂度都为$o(log n)$(恰为线段树上的某个区间),因此维护$FL$的复杂度为$o(log^{2}n)$
对于$ML(l,r,h)$,考虑使用懒标记,即对某个区间打上$h$的懒标记
同时,在打懒标记的同时要计算出此区间的$sum_{i=l}^{r}a_{i}$,使用$FL$即可,复杂度为$o(log n)$(恰为某个区间)
标记的下传不需要按照$FL$,直接暴力下传即可(注意左区间要对$max_{i=mid}^{r-1}b_{i}$取$max$)
最终,打懒标记(以及下传标记,两者其实是相同的)复杂度都是$o(log n)$,一次$ML$操作以及查询的复杂度都是$o(log^{2}n)$
对于$a_{x}<a_{x+1}$的情况做法类似,注意两者的懒标记并不需要合并,直接覆盖即可
综上,即做到了单次$FL(l,r,h)$和$ML(l,r,h)$为$o(log^{2}n)$的复杂度,由于二分总复杂度为$o(qlog^{3}n)$
在精度上,由于比较玄学,这里给出两个或许可行的优化:
1.在查询“最大的$lle x$,使得$b_{l-1}ge a_{x}$”将$a_{x}$减去eps进行查询
2.建议调整eps,可以尝试在$10^{-10}$到$10^{-14}$中调整
第四部分:优化二分
上面的做法在常数足够优秀时可以通过,但事实上还可以做到$o(qlog^{2}n)$
对于$h$的优化比较简单,由于区间已经确定,将这个区间在线段树上划分为$o(log n)$段,假设其从左到右依次为$[l_{1},r_{1}],[l_{2},r_{2}],...,[l_{k},r_{k}]$,并令$mx_{i}=max_{j=l_{i}-1}^{r_{i}-1}b_{j}$
考虑总恰存在$0le ile k$,使得$max_{j=i+1}^{k}mx_{j}le h<max_{j=i}^{k}mx_{j}$(特别的,定义$mx_{0}=infty$)
通过$i$,则有$FL(l,r,h)=FL(l_{i},r_{i},h)+FL(l_{1},r_{i-1},max_{j=i}^{k}mx_{j})+sum_{j=i+1}^{k}(r_{j}-l_{j}+1)h$
(特别的,当$i=0$时前两项都定义为0)
计算$FL(l_{i},r_{i},h)$的复杂度为$o(log n)$,中间的项记作$f_{i-1}$,有$f_{i}=f_{i-1}+FL(l_{i},r_{i},max_{j=i+1}^{k}mx_{j})$,且由于与$h$无关,可以$o(log^{2}n)$预处理出来,最后一项直接后缀和即可
由此,即可$o(log n)$的计算$FL(l,r,h)$,这一部分总复杂度降为$o(log^{2}n)$
(最终的代码就优化到这里,最终开Ofast卡过了)
对于$r$的二分优化,我并没有看懂,下面给出论文中的原文——
我们发现它本质要求任意区间右边接一个线段树节点,询问新区间的$FL$或$FR$函数。我们可以考虑直接新建一个线段树节点存储这个新区间的区间信息,注意每个节点需要存储它的左右儿子。我们发现合并两个区间后询问$FL$和$FR$仍然是$o(log n)$的,于是我们在线段树上二分时维护这个合并出的新节点即可将复杂度降到$o(log^{2}n)$
1 #include<bits/stdc++.h> 2 using namespace std; 3 #define N 100005 4 #define S (N<<2) 5 #define L (k<<1) 6 #define R (L+1) 7 #define mid (l+r>>1) 8 #define ld long double 9 #define eps 1e-14 10 struct Seg{ 11 int k,l,r; 12 }; 13 vector<Seg>v; 14 int n,m,p,x,sum[N]; 15 ld y,Mx,a[N],b[N],mxB[S],sumA[S],type[S],tag[S],mxL[S],mxR[S],FL[S],FR[S],mx[N],f[N]; 16 void up_B(int k){ 17 mxB[k]=max(mxB[L],mxB[R]); 18 } 19 void build_B(int k,int l,int r){ 20 if (l==r){ 21 mxB[k]=b[l]; 22 return; 23 } 24 build_B(L,l,mid); 25 build_B(R,mid+1,r); 26 up_B(k); 27 } 28 void update_B(int k,int l,int r,int x,ld y){ 29 if (l==r){ 30 mxB[k]=y; 31 return; 32 } 33 if (x<=mid)update_B(L,l,mid,x,y); 34 else update_B(R,mid+1,r,x,y); 35 up_B(k); 36 } 37 ld query_B(int k,int l,int r,int x,int y){ 38 if ((l>y)||(x>r))return 0; 39 if ((x<=l)&&(r<=y))return mxB[k]; 40 return max(query_B(L,l,mid,x,y),query_B(R,mid+1,r,x,y)); 41 } 42 int getpre_B(int k,int l,int r,int x,ld y){ 43 if ((l>x)||(mxB[k]<y))return 0; 44 if (l==r)return l; 45 int ans=getpre_B(R,mid+1,r,x,y); 46 if (ans)return ans; 47 return getpre_B(L,l,mid,x,y); 48 } 49 int getnex_B(int k,int l,int r,int x,ld y){ 50 if ((r<x)||(mxB[k]<y))return n; 51 if (l==r)return l; 52 int ans=getnex_B(L,l,mid,x,y); 53 if (ans<n)return ans; 54 return getnex_B(R,mid+1,r,x,y); 55 } 56 ld query_FL(int k,int l,int r,ld x){ 57 if (l==r)return x; 58 if (mxL[k]<=x)return query_FL(L,l,mid,x)+(r-mid)*x; 59 return query_FL(R,mid+1,r,x)+FL[k]; 60 } 61 ld query_FR(int k,int l,int r,ld x){ 62 if (l==r)return x; 63 if (mxR[k]<=x)return query_FR(R,mid+1,r,x)+(mid-l+1)*x; 64 return query_FR(L,l,mid,x)+FR[k]; 65 } 66 ld query_FL(int k,int l,int r,int x,int y){ 67 if ((l>y)||(x>r))return 0; 68 if ((x<=l)&&(r<=y)){ 69 ld ans=query_FL(k,l,r,Mx); 70 Mx=max(Mx,query_B(1,1,n-1,l-1,r-1)); 71 return ans; 72 } 73 ld ans=query_FL(R,mid+1,r,x,y); 74 return ans+query_FL(L,l,mid,x,y); 75 } 76 ld query_FR(int k,int l,int r,int x,int y){ 77 if ((l>y)||(x>r))return 0; 78 if ((x<=l)&&(r<=y)){ 79 ld ans=query_FR(k,l,r,Mx); 80 Mx=max(Mx,query_B(1,1,n-1,l,r)); 81 return ans; 82 } 83 ld ans=query_FR(L,l,mid,x,y); 84 return ans+query_FR(R,mid+1,r,x,y); 85 } 86 ld Query_FL(int x,int y,ld z){ 87 Mx=z; 88 return query_FL(1,1,n,x,y); 89 } 90 ld Query_FR(int x,int y,ld z){ 91 Mx=z; 92 return query_FR(1,1,n,x,y); 93 } 94 void upd_tag(int k,int l,int r,int x,ld y){ 95 type[k]=x,tag[k]=y; 96 if (x==1)sumA[k]=query_FL(k,l,r,y); 97 if (x==2)sumA[k]=query_FR(k,l,r,y); 98 } 99 void up_A(int k){ 100 sumA[k]=sumA[L]+sumA[R]; 101 } 102 void down_A(int k,int l,int r){ 103 if (!type[k])return; 104 if (type[k]==1){ 105 upd_tag(L,l,mid,1,max(tag[k],mxL[k])); 106 upd_tag(R,mid+1,r,1,tag[k]); 107 } 108 if (type[k]==2){ 109 upd_tag(L,l,mid,2,tag[k]); 110 upd_tag(R,mid+1,r,2,max(tag[k],mxR[k])); 111 } 112 type[k]=0; 113 } 114 void upd_mxL(int k,int l,int r){ 115 mxL[k]=query_B(1,1,n-1,mid,r-1); 116 FL[k]=query_FL(L,l,mid,mxL[k]); 117 } 118 void upd_mxR(int k,int l,int r){ 119 mxR[k]=query_B(1,1,n-1,l,mid); 120 FR[k]=query_FR(R,mid+1,r,mxR[k]); 121 } 122 void build_A(int k,int l,int r){ 123 upd_tag(k,l,r,0,0); 124 if (l==r){ 125 sumA[k]=a[l]; 126 return; 127 } 128 build_A(L,l,mid); 129 build_A(R,mid+1,r); 130 upd_mxL(k,l,r),upd_mxR(k,l,r); 131 up_A(k); 132 } 133 void update_A(int k,int l,int r,int x,ld y){ 134 if (l==r)return; 135 down_A(k,l,r); 136 if (x<=mid)update_A(L,l,mid,x,y); 137 else update_A(R,mid+1,r,x,y); 138 upd_mxL(k,l,r),upd_mxR(k,l,r); 139 } 140 ld update_ML(int k,int l,int r,int x,int y,ld z){ 141 if ((l>y)||(x>r))return 0; 142 if ((x<=l)&&(r<=y)){ 143 upd_tag(k,l,r,1,z); 144 return query_B(1,1,n-1,l-1,r-1); 145 } 146 down_A(k,l,r); 147 ld s=update_ML(R,mid+1,r,x,y,z); 148 s=max(s,update_ML(L,l,mid,x,y,max(z,s))); 149 up_A(k); 150 return s; 151 } 152 ld update_MR(int k,int l,int r,int x,int y,ld z){ 153 if ((l>y)||(x>r))return 0; 154 if ((x<=l)&&(r<=y)){ 155 upd_tag(k,l,r,2,z); 156 return query_B(1,1,n-1,l,r); 157 } 158 down_A(k,l,r); 159 ld s=update_MR(L,l,mid,x,y,z); 160 s=max(s,update_MR(R,mid+1,r,x,y,max(z,s))); 161 up_A(k); 162 return s; 163 } 164 ld query_A(int k,int l,int r,int x,int y){ 165 if ((l>y)||(x>r))return 0; 166 if ((x<=l)&&(r<=y))return sumA[k]; 167 down_A(k,l,r); 168 return query_A(L,l,mid,x,y)+query_A(R,mid+1,r,x,y); 169 } 170 ld query(int x){ 171 return query_A(1,1,n,x,x); 172 } 173 void getseg(int k,int l,int r,int x,int y){ 174 if ((l>y)||(x>r))return; 175 if ((x<=l)&&(r<=y)){ 176 v.push_back(Seg{k,l,r}); 177 return; 178 } 179 getseg(L,l,mid,x,y); 180 getseg(R,mid+1,r,x,y); 181 } 182 ld calcL(ld h){ 183 for(int i=(int)v.size()-1;i>=0;i--) 184 if (mx[i]>h)return query_FL(v[i].k,v[i].l,v[i].r,h)+f[i]+sum[i+1]*h; 185 return sum[0]*h; 186 } 187 ld calcR(ld h){ 188 for(int i=1;i<=v.size();i++) 189 if (mx[i]>h)return query_FR(v[i-1].k,v[i-1].l,v[i-1].r,h)+f[i]+sum[i-1]*h; 190 return sum[v.size()]*h; 191 } 192 void update(int x,ld y){ 193 b[x]=y; 194 update_B(1,1,n-1,x,y); 195 update_A(1,1,n,x,y); 196 if ((query(x)<=b[x])&&(query(x+1)<=b[x]))return; 197 if (query(x)>query(x+1)){ 198 int l=getpre_B(1,1,n-1,x-1,query(x)-eps)+1; 199 int i=x+1,j=n; 200 while (i<j){ 201 int k=(i+j+1>>1); 202 if (Query_FL(l,k,query(k))<=query_A(1,1,n,l,k))i=k; 203 else j=k-1; 204 } 205 ld ii=0,jj=100,s=query_A(1,1,n,l,i); 206 v.clear(); 207 getseg(1,1,n,l,i); 208 for(int i=0;i<v.size();i++)mx[i]=query_B(1,1,n-1,v[i].l-1,v[i].r-1); 209 mx[v.size()]=sum[v.size()]=0; 210 for(int i=v.size();i;i--){ 211 mx[i-1]=max(mx[i],mx[i-1]); 212 sum[i-1]=sum[i]+v[i-1].r-v[i-1].l+1; 213 } 214 f[0]=0; 215 for(int i=0;i<v.size();i++)f[i+1]=f[i]+query_FL(v[i].k,v[i].l,v[i].r,mx[i+1]); 216 while (jj-ii>eps){ 217 ld k=(ii+jj)/2; 218 if (calcL(k)<=s)ii=k; 219 else jj=k; 220 } 221 update_ML(1,1,n,l,i,(ii+jj)/2); 222 } 223 else{ 224 int r=getnex_B(1,1,n-1,x+1,query(x+1)-eps); 225 int i=1,j=x; 226 while (i<j){ 227 int k=(i+j>>1); 228 if (Query_FR(k,r,query(k))<=query_A(1,1,n,k,r))j=k; 229 else i=k+1; 230 } 231 ld ii=0,jj=100,s=query_A(1,1,n,i,r); 232 v.clear(); 233 getseg(1,1,n,i,r); 234 mx[0]=sum[0]=0; 235 for(int i=0;i<v.size();i++){ 236 mx[i+1]=max(mx[i],query_B(1,1,n-1,v[i].l,v[i].r)); 237 sum[i+1]=sum[i]+v[i].r-v[i].l+1; 238 } 239 f[v.size()]=0; 240 for(int i=(int)v.size()-1;i>=0;i--)f[i]=f[i+1]+query_FR(v[i].k,v[i].l,v[i].r,mx[i]); 241 while (jj-ii>eps){ 242 ld k=(ii+jj)/2; 243 if (calcR(k)<=s)ii=k; 244 else jj=k; 245 } 246 update_MR(1,1,n,i,r,(ii+jj)/2); 247 } 248 } 249 int main(){ 250 scanf("%d%d",&n,&m); 251 for(int i=1;i<=n;i++)scanf("%Lf",&a[i]); 252 for(int i=1;i<n;i++)scanf("%Lf",&b[i]); 253 build_B(1,1,n-1); 254 build_A(1,1,n); 255 for(int i=1;i<=m;i++){ 256 scanf("%d%d",&p,&x); 257 if (p==1){ 258 scanf("%Lf",&y); 259 if (b[x]>y)update(x,y); 260 } 261 if (p==2)printf("%.8Lf ",query(x)); 262 } 263 for(int i=1;i<=n;i++)printf("%.8Lf ",query(i)); 264 }