题解：B3646 数列前缀和 3

· · 题解

B3646 数列前缀和 3

背景

不会矩阵求逆，遂选择线段树做法，做完发现 TLE 了，经过疯狂的常数优化后，通过了这道题。

前置知识

矩阵乘法，线段树。

分析

这道题目给了我们一些矩阵，让我们多次求出一个给定区间内所有矩阵的乘积。

这道题当然可以直接用矩阵前缀积来做，但是因为要取模，所以需要求出矩阵逆，这样问题就变得复杂了。怎么办呢？线段树！

我们在线段树的每个叶子节点上面存一个矩阵，其余的存储区间矩阵积。最开始把所有的矩阵存进去，然后每次查询的时候按照题目要求异或即可。

需要注意的是，线段树查询操作累乘答案的时候，答案的初值应为单位元矩阵，否则会影响答案。

还需要一个常数优化的技巧：循环展开。

循环展开，英文中称 （Loop unwinding 或 loop unrolling） ，是一种牺牲程序的尺寸来加快程序的执行速度的优化方法。

循环展开可以省去 for 循环每次的判断与修改操作，对于多次调用循环的程序，把循环展开可以有效卡常。

其余的就是一些常见的小细节了，放代码。

Code

#include<bits/stdc++.h>
using namespace std;
#define ll unsigned long long
#define ls (p<<1)
#define rs (p<<1|1)
inline int read()//快读 
{
    int s=0;char c=getchar();
    for(;!isdigit(c);c=getchar());
    for(;isdigit(c);c=getchar())
     s=(s<<1)+(s<<3)+(c^48);
    return s;
}
int const N=1e6+9;
int const mod=1145141;
int n,k,q,ans;
ll tmp[3][3];
struct matrix{
    ll a[3][3];
    inline void operator *=(const matrix &B)//矩阵乘法 
    {
        memcpy(tmp,a,sizeof a);
        if(k==3)//k=3时的循环展开 
        {
            a[0][0]=(tmp[0][0]*B.a[0][0]+tmp[0][1]*B.a[1][0]+tmp[0][2]*B.a[2][0])%mod;
            a[0][1]=(tmp[0][0]*B.a[0][1]+tmp[0][1]*B.a[1][1]+tmp[0][2]*B.a[2][1])%mod;
            a[0][2]=(tmp[0][0]*B.a[0][2]+tmp[0][1]*B.a[1][2]+tmp[0][2]*B.a[2][2])%mod;
            a[1][0]=(tmp[1][0]*B.a[0][0]+tmp[1][1]*B.a[1][0]+tmp[1][2]*B.a[2][0])%mod;
            a[1][1]=(tmp[1][0]*B.a[0][1]+tmp[1][1]*B.a[1][1]+tmp[1][2]*B.a[2][1])%mod;
            a[1][2]=(tmp[1][0]*B.a[0][2]+tmp[1][1]*B.a[1][2]+tmp[1][2]*B.a[2][2])%mod;
            a[2][0]=(tmp[2][0]*B.a[0][0]+tmp[2][1]*B.a[1][0]+tmp[2][2]*B.a[2][0])%mod;
            a[2][1]=(tmp[2][0]*B.a[0][1]+tmp[2][1]*B.a[1][1]+tmp[2][2]*B.a[2][1])%mod;
            a[2][2]=(tmp[2][0]*B.a[0][2]+tmp[2][1]*B.a[1][2]+tmp[2][2]*B.a[2][2])%mod;
        }else {//k=2 
            a[0][0]=(tmp[0][0]*B.a[0][0]+tmp[0][1]*B.a[1][0])%mod;
            a[0][1]=(tmp[0][0]*B.a[0][1]+tmp[0][1]*B.a[1][1])%mod;
            a[1][0]=(tmp[1][0]*B.a[0][0]+tmp[1][1]*B.a[1][0])%mod;
            a[1][1]=(tmp[1][0]*B.a[0][1]+tmp[1][1]*B.a[1][1])%mod;
        }
    }
    inline int xor_()//矩阵异或 
    {
        int res=0;
        if(k==3)//还是循环展开 
        {
            res^=a[0][0];
            res^=a[0][1];
            res^=a[0][2];

            res^=a[1][0];
            res^=a[1][1];
            res^=a[1][2];

            res^=a[2][0];
            res^=a[2][1];
            res^=a[2][2];
        }else{
            res^=a[0][0];
            res^=a[0][1];

            res^=a[1][0];
            res^=a[1][1];
        }
        return res;
    }
}tr[N<<2],a[N];
void build(int l,int r,int p)//建树 
{
    if(l==r)
    {
        tr[p]=a[l];
        return;
    }
    int mid=(l+r)>>1;
    build(l,mid,ls);
    build(mid+1,r,rs);
    tr[p]=tr[ls];
    tr[p]*=tr[rs];
}
matrix query(int nl,int nr,int l,int r,int p)//查询 
{
    if(nl<=l&&r<=nr)
     return tr[p];
    int mid=(l+r)>>1;
    matrix res;memset(res.a,0,sizeof(res.a));
    res.a[0][0]=res.a[1][1]=res.a[2][2]=1;
    if(nl<=mid)res*=query(nl,nr,l,mid,ls);
    if(nr>mid)res*=query(nl,nr,mid+1,r,rs);
    return res;
}
signed main()
{
    n=read();k=read();q=read();
    for(int t=1;t<=n;t++)
    {
        if(k==3)
        {
            a[t].a[0][0]=read();
            a[t].a[0][1]=read();
            a[t].a[0][2]=read();

            a[t].a[1][0]=read();
            a[t].a[1][1]=read();
            a[t].a[1][2]=read();

            a[t].a[2][0]=read();
            a[t].a[2][1]=read();
            a[t].a[2][2]=read();
        }else{
            a[t].a[0][0]=read();
            a[t].a[0][1]=read();

            a[t].a[1][0]=read();
            a[t].a[1][1]=read();
        }
    }
    build(1,n,1);
    int x,y;
    matrix z;
    while(q--)
    {
        x=read();y=read();
        z=query(x,y,1,n,1);
        ans^=z.xor_();
    }
    printf("%d",ans);
}

特别鸣谢

@lkwbian