「LibreOJ NOI Round #2」单枪匹马

题目

通过这道题成功发现我不会矩乘

答案是一个连分数,看起来不像是一般的数据结构能做的样子,设(x_{l,r},y_{l,r})分别表示([l,r])询问的分子和分母

于是有

[frac{x_{l,r}}{y_{l,r}}=a_{l}+frac{y_{l+1,r}}{x_{l+1,r}}=frac{y_{l+1,r}+a_l imes x_{l+1,r}}{x_{l+1,r}} ]

这么鬼畜的东西就应该写成矩阵的形式

于是

[egin{bmatrix} a_l& 1\ 1&0 end{bmatrix} imes egin{bmatrix} x_{l+1,r}\ y_{l+1,r} end{bmatrix}=egin{bmatrix} y_{l+1,r}+a_l imes x_{l+1,r}\ x_{l+1,r} end{bmatrix}= egin{bmatrix} x_{l,r}\ y_{l,r} end{bmatrix} ]

这是一个(2 imes 2)的矩阵乘一个(2 imes 1)的向量,所以必须是矩阵在前去乘向量其实我之前一直把向量写前面

于是对于(l,r)这样的一组询问我们要求的就是(egin{bmatrix} a_r\ 1 end{bmatrix})这样一个向量,去乘(r-1)(l)这些矩阵,而且必须是从(r-1)乘到(l),由于乘完一个矩阵得到的还是一个向量,所以向量还是放在后面,于是根据结合律前面就是从(l)乘到(r-1)

就是

[egin{bmatrix} a_l& 1\ 1&0 end{bmatrix} imes egin{bmatrix} a_{l+1}& 1\ 1&0 end{bmatrix} imes egin{bmatrix} a_{l+2}& 1\ 1&0 end{bmatrix} imes ... imes egin{bmatrix} a_{r-1}& 1\ 1&0 end{bmatrix} imes egin{bmatrix} a_r\ 1 end{bmatrix} ]

于是我们需要查询一个区间内矩阵的乘积,数据范围(10^6),看起来不太好用带(log)的东西维护,于是考虑强行求逆

(2 imes 2)的逆矩阵手推一下就推出来了,就是(egin{bmatrix} 0& 1\ 1&-a_{l} end{bmatrix}),于是我们再维护一个逆矩阵的前缀积就可以快速查询一个区间了

但是逆矩阵的前缀积需要维护从(i)乘到(1)的积,矩阵的前缀积需要维护从(1)(i)的积,查询的时候也是那逆矩阵的前缀积乘上矩阵前缀积,这样前面才能两两消掉

代码

#include<bits/stdc++.h>
#define re register
#define LL long long
#define max(a,b) ((a)>(b)?(a):(b))
#define min(a,b) ((a)<(b)?(a):(b))
inline int read() {
	char c=getchar();int x=0;while(c<'0'||c>'9') c=getchar();
	while(c>='0'&&c<='9') x=(x<<3)+(x<<1)+c-48,c=getchar();return x;
}
const int maxn=1e6+5;
const int mod=998244353;
struct mat{int a[2][2];};
inline int qm(int x) {return x>=mod?x-mod:x;}
inline mat operator*(mat a,mat b) {
	mat c;
	c.a[0][0]=qm(1ll*a.a[0][0]*b.a[0][0]%mod+1ll*a.a[0][1]*b.a[1][0]%mod);
        c.a[0][1]=qm(1ll*a.a[0][0]*b.a[0][1]%mod+1ll*a.a[0][1]*b.a[1][1]%mod);
        c.a[1][0]=qm(1ll*a.a[1][0]*b.a[0][0]%mod+1ll*a.a[1][1]*b.a[1][0]%mod);
        c.a[1][1]=qm(1ll*a.a[1][0]*b.a[0][1]%mod+1ll*a.a[1][1]*b.a[1][1]%mod);
        return c;
} 
inline void out(mat c) {
	printf("%d %d
",c.a[0][0],c.a[0][1]);
	printf("%d %d
",c.a[1][0],c.a[1][1]);
	puts("");
}
mat pre[maxn],inv[maxn];
int n,m,type,a[maxn];
int main() {
	n=read(),m=read(),type=read();
	pre[0].a[0][0]=pre[0].a[1][1]=inv[0].a[0][0]=inv[0].a[1][1]=1;
	for(re int i=1;i<=n;i++) pre[i].a[0][1]=pre[i].a[1][0]=1;
	for(re int i=1;i<=n;i++) inv[i].a[0][1]=inv[i].a[1][0]=1;
	for(re int i=1;i<=n;i++) pre[i].a[0][0]=a[i]=read(),inv[i].a[1][1]=qm(mod-a[i]);
	for(re int i=1;i<=n;i++) pre[i]=pre[i-1]*pre[i],inv[i]=inv[i]*inv[i-1];
	int lstx=0,lsty=0,op,x,l,r;
	while(m--) {
		op=read();
		if(op==1) {
			x=read();
			if(type) x^=(lstx^lsty);a[++n]=x;
			pre[n].a[0][1]=pre[n].a[1][0]=inv[n].a[0][1]=inv[n].a[1][0]=1;
			pre[n].a[0][0]=x,inv[n].a[1][1]=qm(mod-x);
			pre[n]=pre[n-1]*pre[n],inv[n]=inv[n]*inv[n-1];
		}
		if(op==2) {
			l=read(),r=read();
			if(type) l^=(lstx^lsty),r^=(lstx^lsty);
			if(l==r) lstx=a[l],lsty=1;
			else {
				mat c=inv[l-1]*pre[r-1];
				lstx=qm(1ll*a[r]*c.a[0][0]%mod+c.a[0][1]);
				lsty=qm(1ll*a[r]*c.a[1][0]%mod+c.a[1][1]);
			}
			printf("%d %d
",lstx,lsty);
		}
	}
	return 0;
}

原文地址:https://www.cnblogs.com/asuldb/p/11542312.html