题解：[NOIP2022] 比赛

题目传送门

题意分析

原题面等价于给定序列 \(a_1,a_2,a_3,\cdots,a_n\) 和 \(b_1,b_2,b_3,\cdots,b_n\) 和 \(q\) 次询问，每次询问给定 \(l,r\)，求：

\[\sum_{p=l}^r\sum_{q=p}^r\max_{i=p}^qa_i\max_{i=p}^qb_i\bmod2^{64} \]

ST 表暴力做是 \(\mathcal O\left(qn^2\right)\) 的，可以得到 8pts 的高分。

设：

\[\begin{aligned} x_{l,r}&=\max_{i=l}^ra_i\\ y_{l,r}&=\max_{i=l}^rb_i \end{aligned} \]

则，答案就是：

\[\sum_{p=l}^r\sum_{q=l}^rx_{p,q}y_{p,q} \]

考虑对 \(r\) 进行扫描线，可以得到转移：

\[\begin{aligned} x_{l,r}&=\max(x_{l,r-1},a_r)\\ y_{l,r}&=\max(y_{l,r-1},b_r)\\ \end{aligned} \]

对于询问 \((l_0,r_0)\)，答案为：

\[\sum_{r=l_0}^{r_0}\sum_{l=l_0}^rx_{l,r}y_{l,r} \]

因此我们在扫描 \(r\) 的过程中，需要维护 \(\displaystyle\sum_{l=l_0}^rx_{l,r}y_{l,r}\) 的历史和。

为了便于讨论，钦定 \(r\) 一定，记为 \(x_l,y_l\)，修改就直接做 \(\operatorname{chkmax}\)。钦定 \(r<l\) 时 \(x_l=y_l=0\)，我们要求的就是 \(\displaystyle\sum_{i=l_0}^{r_0}x_iy_i\) 的 \([l_0,r_0]\) 版本的历史版本和。

离线后可以差分处理，但是因为我们钦定了初始为 \(0\)，因此差不差分都一样。扫到 \(l_0\) 之前 \(l_0\) 右边都是 \(0\)，直接求所有版本的历史和即可。

暴力做可以得到 20pts，考虑数据结构维护，区间操作容易想到线段树。

---

如果是线段树维护 \(\operatorname{chkmax}\) 的话是愚蠢的，因为这样要写吉司机线段树，很难写还有维护很多信息。

因为 \(a,b\) 已知，因此可以提前预处理每一个 \(a_i,b_i\) 能够影响的区间，从而将线段树操作转化为区间覆盖。具体而言，以 \(a_i\) 为例，其覆盖范围为 \([j+1,i]\)，\(j\) 为满足 \(j<i\) 的最后一个大于 \(a_i\) 的数。换而言之，实际上 \(x,y\) 分别是 \(a,b\) 的动态后缀最大值。这很容易用单调栈得到。

之后我们需要的线段树操作有：

• \(x_i,y_i\) 区间覆盖
• 查询 \(\displaystyle\sum x_iy_i\) 历史和

重点显然在于 \(\displaystyle\sum x_iy_i\)。考虑线段树维护历史和。

历史和线段树

先考虑矩阵乘法实现的历史和线段树。

矩阵乘法

记 \(\textit{hxy}\) 为 \(\displaystyle\sum x_iy_i\) 历史和，简记 \(\textit{xy}\) 表示 \(x_iy_i\)，\(x\) 为 \(x_i\)，\(y\) 为 \(y_i\)。因为涉及到区间赋值，需要用到区间长度判断新的区间和，因此把区间长度 \(d\) 丢入矩阵。但是这只是为了转移方便，事实上对于所有的线段树节点，其区间长度 \(d\) 都为定值。

\[\begin{bmatrix} \textit{hxy}\\ \displaystyle\sum xy\\ \displaystyle\sum x\\ \displaystyle\sum y\\ d \end{bmatrix} \]

对于求历史和，有转移矩阵：

\[\begin{bmatrix} 1&1&0&0&0\\ 0&1&0&0&0\\ 0&0&1&0&0\\ 0&0&0&1&0\\ 0&0&0&0&1 \end{bmatrix} \begin{bmatrix} \textit{hxy}\\ \displaystyle\sum xy\\ \displaystyle\sum x\\ \displaystyle\sum y\\ d \end{bmatrix} = \begin{bmatrix} \displaystyle \textit{hxy}+\sum xy\\ \displaystyle\sum xy\\ \displaystyle\sum x\\ \displaystyle\sum y\\ d \end{bmatrix} \]

之后考虑 \(x\) 区间覆盖为 \(k\)，有：

\[\begin{bmatrix} 1&0&0&0&0\\ 0&0&0&k&0\\ 0&0&0&0&k\\ 0&0&0&1&0\\ 0&0&0&0&1 \end{bmatrix} \begin{bmatrix} \textit{hxy}\\ \displaystyle\sum xy\\ \displaystyle\sum x\\ \displaystyle\sum y\\ d \end{bmatrix} = \begin{bmatrix} \displaystyle \textit{hxy}\\ \displaystyle k\sum y\\ dk\\ \displaystyle\sum y\\ d \end{bmatrix} \]

同理还有 \(y\) 区间覆盖为 \(k\)，有：

\[\begin{bmatrix} 1&0&0&0&0\\ 0&0&k&0&0\\ 0&0&1&0&0\\ 0&0&0&0&k\\ 0&0&0&0&1 \end{bmatrix} \begin{bmatrix} \textit{hxy}\\ \displaystyle\sum xy\\ \displaystyle\sum x\\ \displaystyle\sum y\\ d \end{bmatrix} = \begin{bmatrix} \displaystyle \textit{hxy}\\ \displaystyle k\sum x\\ \displaystyle\sum x\\ dk\\ d \end{bmatrix} \]

矩阵乘法不用动脑子的代价就是常数大吗。——xyh

这样做时间复杂度 \(\mathcal O(q\log n)\)，但是矩阵乘法会有 \(5^3\) 的巨大常数，所以下面这份代码得到了 36pts 的高分：

朴素矩阵乘法参考代码

//#include<bits/stdc++.h>
#include<algorithm>
#include<iostream>
#include<cstring>
#include<iomanip>
#include<cstdio>
#include<string>
#include<vector>
#include<cmath>
#include<ctime>
#include<deque>
#include<queue>
#include<stack>
#include<list>
using namespace std;
constexpr const int N=2.5e5,Q=2.5e5;
typedef unsigned long long ull;
int n,q,a[N+1],b[N+1],preA[N+1],preB[N+1];
ull ans[N+1]; 
vector<pair<int,int>>f[N+1];
void pre(){vector<int>s;for(int i=1;i<=n;i++){while(s.size()&&a[s.back()]<=a[i]){s.pop_back();}if(s.size()){preA[i]=s.back();}s.push_back(i);}s.resize(0);for(int i=1;i<=n;i++){while(s.size()&&b[s.back()]<=b[i]){s.pop_back();}if(s.size()){preB[i]=s.back();}s.push_back(i);}
}
struct Matrix{int n,m;ull a[5][5];Matrix(int nn=-1,int mm=-1){if(mm==-1){mm=nn;}if(nn!=-1){n=nn,m=mm;}else{n=0,m=0;}memset(a,0,sizeof(a));}Matrix(initializer_list<initializer_list<ull>> x){n=x.size();m=0;for(auto &i:x){m=max(m,(int)i.size());}int pi=0;for(auto &i:x){int pj=0;for(auto &j:i){a[pi][pj]=j;pj++;}pi++;}}ull* operator [](const int &x){return a[x];}void unit(){memset(a,0,sizeof(a));for(int i=0;i<n;i++){a[i][i]=1;}} void print(){for(int i=0;i<n;i++){for(int j=0;j<m;j++){cout<<a[i][j]<<' ';}cout<<endl;}}
};
Matrix operator *(Matrix A,Matrix B){Matrix C(A.n,B.m);for(int i=0;i<A.n;i++){for(int j=0;j<B.m;j++){for(int k=0;k<A.m;k++){C[i][j]+=A[i][k]*B[k][j];}}}return C;
}
Matrix& operator *=(Matrix &A,Matrix B){return A=A*B;
}
Matrix operator +(Matrix A,Matrix B){Matrix C(A.n,A.m);for(int i=0;i<A.n;i++){for(int j=0;j<A.m;j++){C[i][j]=A[i][j]+B[i][j];} } return C;
}
Matrix& operator +=(Matrix &A,Matrix B){for(int i=0;i<A.n;i++){for(int j=0;j<A.m;j++){A[i][j]+=B[i][j];}}return A;
}
struct segTree{struct node{int l,r;Matrix value,tag;}t[N<<2|1];void up(int p){t[p].value=t[p<<1].value+t[p<<1|1].value;}void build(int p,int l,int r){t[p].l=l,t[p].r=r;t[p].value=Matrix(5,1);t[p].tag=Matrix(5);t[p].tag.unit();if(l==r){t[p].value[4][0]=1;return;}int mid=l+r>>1;build(p<<1,l,mid);build(p<<1|1,mid+1,r);up(p);}void down(int p){t[p<<1].value=t[p].tag*t[p<<1].value;t[p<<1].tag=t[p].tag*t[p<<1].tag;t[p<<1|1].value=t[p].tag*t[p<<1|1].value;t[p<<1|1].tag=t[p].tag*t[p<<1|1].tag;t[p].tag.unit(); }void setX(int p,int l,int r,int k){if(l<=t[p].l&&t[p].r<=r){t[p].value=Matrix({{1,0,0,0,0},{0,0,0,k,0},{0,0,0,0,k},{0,0,0,1,0},{0,0,0,0,1}})*t[p].value;t[p].tag=Matrix({{1,0,0,0,0},{0,0,0,k,0},{0,0,0,0,k},{0,0,0,1,0},{0,0,0,0,1}})*t[p].tag;return;}down(p);if(l<=t[p<<1].r){setX(p<<1,l,r,k);}if(t[p<<1|1].l<=r){setX(p<<1|1,l,r,k);}up(p);}void setY(int p,int l,int r,int k){if(l<=t[p].l&&t[p].r<=r){t[p].value=Matrix({{1,0,0,0,0},{0,0,k,0,0},{0,0,1,0,0},{0,0,0,0,k},{0,0,0,0,1}})*t[p].value;t[p].tag=Matrix({{1,0,0,0,0},{0,0,k,0,0},{0,0,1,0,0},{0,0,0,0,k},{0,0,0,0,1}})*t[p].tag;return;}down(p);if(l<=t[p<<1].r){setY(p<<1,l,r,k);}if(t[p<<1|1].l<=r){setY(p<<1|1,l,r,k);}up(p);}void count(){t[1].value=Matrix({{1,1,0,0,0},{0,1,0,0,0},{0,0,1,0,0},{0,0,0,1,0},{0,0,0,0,1}})*t[1].value;t[1].tag=Matrix({{1,1,0,0,0},{0,1,0,0,0},{0,0,1,0,0},{0,0,0,1,0},{0,0,0,0,1}})*t[1].tag;}ull query(int p,int l,int r){if(l<=t[p].l&&t[p].r<=r){return t[p].value[0][0];}down(p);ull ans=0;if(l<=t[p<<1].r){ans+=query(p<<1,l,r);}if(t[p<<1|1].l<=r){ans+=query(p<<1|1,l,r);}return ans;}
}t;
int main(){/*freopen("test.in","r",stdin);freopen("test.out","w",stdout);*/ios::sync_with_stdio(false);cin.tie(0);cout.tie(0);int T;cin>>T>>n;for(int i=1;i<=n;i++){cin>>a[i];}for(int i=1;i<=n;i++){cin>>b[i];}cin>>q;for(int i=1;i<=q;i++){int l,r;cin>>l>>r;f[r].push_back({l,i});}pre();t.build(1,1,n);for(int r=1;r<=n;r++){t.setX(1,preA[r]+1,r,a[r]);t.setY(1,preB[r]+1,r,b[r]);t.count();for(auto [l,id]:f[r]){ans[id]=t.query(1,l,r);}}for(int i=1;i<=q;i++){cout<<ans[i]<<'\n';}cout.flush();/*fclose(stdin);fclose(stdout);*/return 0;
}

拆矩阵乘法

那么，如何进一步优化呢？

考虑矩阵乘法维护转移，实际上就是 \(\mathcal O\left(n^2\right)\) 枚举了每一个标记之间的线性转移关系。

但是有一些标记之间毫无关联，例如本题中的 \(d\) 与其余四个标记毫无关系；那么我们再在矩阵乘法里维护它们的转移，就会浪费大量时间，导致 TLE。换而言之，我们可以确定矩阵中的一些项，使其永远为 \(0\)。

因此，从标记间依赖关系来拆矩阵乘法即可。

可以考虑列表，用「\」表示没有依赖关系。下表表示用纵列的标记得到横行的标记的关系。例如第一列第二行表示 \(\textit{hxy}\) 不能更新 \(\displaystyle\sum xy\)，即得到 \(\displaystyle\sum xy\) 的转移和 \(\textit{hxy}\) 无关。（其实就是矩阵乘法的顺序。）

观察矩阵转移可得：

	\(\textit{hxy}\)	\(\displaystyle\sum xy\)	\(\displaystyle\sum x\)	\(\displaystyle\sum y\)	\(d\)
\(\textit{hxy}\)
\(\displaystyle\sum xy\)	\
\(\displaystyle\sum x\)	\	\		\
\(\displaystyle\sum y\)	\	\	\
\(d\)	\	\	\	\

需要注意依赖关系的依赖。例如 \(d\) 能更新 \(\displaystyle\sum x,\sum y\)，最终更新 \(\displaystyle\sum xy\)。

因此等价于可以确定转移矩阵中的一些 \(0\)：

\[\begin{bmatrix} &&&&\\ 0&&&&\\ 0&0&&0&\\ 0&0&0&&\\ 0&0&0&0&\ \end{bmatrix} \]

一种更聪明的做法

其实可以发现所有的转移矩阵都是主对角线为分割的上三角矩阵。

因此我们仅仅维护上三角即可。

因此可以手动这些非 \(0\) 的维护转移标记：

struct Tag{ull a11,a12,a13,a14,a15,a22,a23,a24,a25,a33,a35,a44,a45,a55;void unit(){a11=a22=a33=a44=a55=1;a12=a13=a14=a15=a23=a24=a25=a35=a45=0;}
};
Tag operator *(Tag a,Tag b){Tag c;c.a11=a.a11*b.a11;c.a12=a.a11*b.a12+a.a12*b.a22;c.a13=a.a11*b.a13+a.a12*b.a23+a.a13*b.a33;c.a14=a.a11*b.a14+a.a12*b.a24+a.a14*b.a44;c.a15=a.a11*b.a15+a.a12*b.a25+a.a13*b.a35+a.a14*b.a45+a.a15*b.a55;c.a22=a.a22*b.a22;c.a23=a.a22*b.a23+a.a23*b.a33;c.a24=a.a22*b.a24+a.a24*b.a44;c.a25=a.a22*b.a25+a.a23*b.a35+a.a24*b.a45+a.a25*b.a55;c.a33=a.a33*b.a33;c.a35=a.a33*b.a35+a.a35*b.a55;c.a44=a.a44*b.a44;c.a45=a.a44*b.a45+a.a45*b.a55;c.a55=a.a55*b.a55;return c;
}

以及线段树节点信息矩阵（向量）：

struct Value{ull a1,a2,a3,a4,a5;
};
Value operator *(Tag a,Value b){Value c;c.a1=a.a11*b.a1+a.a12*b.a2+a.a13*b.a3+a.a14*b.a4+a.a15*b.a5;c.a2=a.a22*b.a2+a.a23*b.a3+a.a24*b.a4+a.a25*b.a5;c.a3=a.a33*b.a3+a.a35*b.a5;c.a4=a.a44*b.a4+a.a45*b.a5;c.a5=a.a55*b.a5;return c;
}

时间复杂度

拆完矩阵后，常数从巨大的 \(5^3\) 降低到了 \(14\)，可以接受。

实际上还可以将本质不同系数个数从 \(14\) 继续缩减，但是没必要。首先可以发现这里面有一些恒为 \(1\)。

其次，\(d\rightarrow\displaystyle\sum x\) 的贡献其实和 \(\displaystyle\sum y\rightarrow\sum xy\) 的贡献相同。因此可以去除一些，但是没有必要。

AC 代码

//#include<bits/stdc++.h>
#include<algorithm>
#include<iostream>
#include<cstring>
#include<iomanip>
#include<cstdio>
#include<string>
#include<vector>
#include<cmath>
#include<ctime>
#include<deque>
#include<queue>
#include<stack>
#include<list>
using namespace std;
constexpr const int N=2.5e5,Q=2.5e5;
typedef unsigned long long ull;
int n,q,a[N+1],b[N+1],preA[N+1],preB[N+1];
ull ans[N+1]; 
vector<pair<int,int>>f[N+1];
void pre(){vector<int>s;for(int i=1;i<=n;i++){while(s.size()&&a[s.back()]<=a[i]){s.pop_back();}if(s.size()){preA[i]=s.back();}s.push_back(i);}s.resize(0);for(int i=1;i<=n;i++){while(s.size()&&b[s.back()]<=b[i]){s.pop_back();}if(s.size()){preB[i]=s.back();}s.push_back(i);}
}
struct Value{ull a1,a2,a3,a4,a5;
};
struct Tag{ull a11,a12,a13,a14,a15,a22,a23,a24,a25,a33,a35,a44,a45,a55;void unit(){a11=a22=a33=a44=a55=1;a12=a13=a14=a15=a23=a24=a25=a35=a45=0;}
};
Tag operator *(Tag a,Tag b){Tag c;c.a11=a.a11*b.a11;c.a12=a.a11*b.a12+a.a12*b.a22;c.a13=a.a11*b.a13+a.a12*b.a23+a.a13*b.a33;c.a14=a.a11*b.a14+a.a12*b.a24+a.a14*b.a44;c.a15=a.a11*b.a15+a.a12*b.a25+a.a13*b.a35+a.a14*b.a45+a.a15*b.a55;c.a22=a.a22*b.a22;c.a23=a.a22*b.a23+a.a23*b.a33;c.a24=a.a22*b.a24+a.a24*b.a44;c.a25=a.a22*b.a25+a.a23*b.a35+a.a24*b.a45+a.a25*b.a55;c.a33=a.a33*b.a33;c.a35=a.a33*b.a35+a.a35*b.a55;c.a44=a.a44*b.a44;c.a45=a.a44*b.a45+a.a45*b.a55;c.a55=a.a55*b.a55;return c;
}
Value operator *(Tag a,Value b){Value c;c.a1=a.a11*b.a1+a.a12*b.a2+a.a13*b.a3+a.a14*b.a4+a.a15*b.a5;c.a2=a.a22*b.a2+a.a23*b.a3+a.a24*b.a4+a.a25*b.a5;c.a3=a.a33*b.a3+a.a35*b.a5;c.a4=a.a44*b.a4+a.a45*b.a5;c.a5=a.a55*b.a5;return c;
}
struct segTree{struct node{int l,r;Value value;Tag tag;}t[N<<2|1];void up(int p){t[p].value.a1=t[p<<1].value.a1+t[p<<1|1].value.a1;t[p].value.a2=t[p<<1].value.a2+t[p<<1|1].value.a2;t[p].value.a3=t[p<<1].value.a3+t[p<<1|1].value.a3;t[p].value.a4=t[p<<1].value.a4+t[p<<1|1].value.a4;t[p].value.a5=t[p<<1].value.a5+t[p<<1|1].value.a5;}void build(int p,int l,int r){t[p].l=l,t[p].r=r;t[p].tag.unit();if(l==r){t[p].value.a5=1;return;}int mid=l+r>>1;build(p<<1,l,mid);build(p<<1|1,mid+1,r);up(p);}void down(int p){t[p<<1].value=t[p].tag*t[p<<1].value;t[p<<1].tag=t[p].tag*t[p<<1].tag;t[p<<1|1].value=t[p].tag*t[p<<1|1].value;t[p<<1|1].tag=t[p].tag*t[p<<1|1].tag;t[p].tag.unit(); }void setX(int p,int l,int r,int k){if(l<=t[p].l&&t[p].r<=r){Tag pl;pl.unit();pl.a22=pl.a33=0;pl.a24=pl.a35=k;t[p].value=pl*t[p].value;t[p].tag=pl*t[p].tag;return;}down(p);if(l<=t[p<<1].r){setX(p<<1,l,r,k);}if(t[p<<1|1].l<=r){setX(p<<1|1,l,r,k);}up(p);}void setY(int p,int l,int r,int k){if(l<=t[p].l&&t[p].r<=r){Tag pl;pl.unit();pl.a22=pl.a44=0;pl.a23=pl.a45=k;t[p].value=pl*t[p].value;t[p].tag=pl*t[p].tag;return;}down(p);if(l<=t[p<<1].r){setY(p<<1,l,r,k);}if(t[p<<1|1].l<=r){setY(p<<1|1,l,r,k);}up(p);}void count(){Tag pl;pl.unit();pl.a12=1;t[1].value=pl*t[1].value;t[1].tag=pl*t[1].tag;}ull query(int p,int l,int r){if(l<=t[p].l&&t[p].r<=r){return t[p].value.a1;}down(p);ull ans=0;if(l<=t[p<<1].r){ans+=query(p<<1,l,r);}if(t[p<<1|1].l<=r){ans+=query(p<<1|1,l,r);}return ans;}
}t;
main(){/*freopen("test.in","r",stdin);freopen("test.out","w",stdout);*/ios::sync_with_stdio(false);cin.tie(0);cout.tie(0);int T;cin>>T>>n;for(int i=1;i<=n;i++){cin>>a[i];}for(int i=1;i<=n;i++){cin>>b[i];}cin>>q;for(int i=1;i<=q;i++){int l,r;cin>>l>>r;f[r].push_back({l,i});}pre();t.build(1,1,n);for(int r=1;r<=n;r++){t.setX(1,preA[r]+1,r,a[r]);t.setY(1,preB[r]+1,r,b[r]);t.count();for(auto [l,id]:f[r]){ans[id]=t.query(1,l,r);}}for(int i=1;i<=q;i++){cout<<ans[i]<<'\n';}cout.flush();/*fclose(stdin);fclose(stdout);*/return 0;
}