河南省多校联盟二-F 线段树+矩阵

---恢复内容开始---

1284: SP教数学

时间限制: 2 秒  内存限制: 128 MB
提交: 24  解决: 4

题目描述

输入

输出

对于每组数据的2操作，输出一行对1e9 + 7取模的答案

样例输入

7 4
2 2 1 1 3 3 2
2 1 5
2 6 7
1 3 4 3
2 6 6

样例输出

6
3
2

提示

1 <=  n ,m <=10^5

一道很有趣的ST的题目，有趣在维护节点时用到了矩阵运算，减少了计算量。

首先对于矩阵运算法则，百度百科:

基本性质

1. 乘法结合律： (AB)C=A(BC)．^[2] 
2. 乘法左分配律：(A+B)C=AC+BC^[2] 
3. 乘法右分配律：C(A+B)=CA+CB^[2] 
4. 对数乘的结合性k(AB）=(kA)B=A(kB）．
5. 转置 (AB)^T=B^TA^T．
6. 矩阵乘法一般不满足交换律^[3] 
   。

题目的重点是如何快速求出∑^r_i=l f(i) , 其中f(i)表示第i个斐波那契数，(f1=f2=1)

这里的i会更新，+x操作，那么如何快速的计算这个值呢，一项一项利用快速幂? len*log（n)的复杂度难免有些太高了。

我们可以利用当前已知的区间和和增量一次性计算得到新的区间和，这里用到一些公式。

假设当前区间有三个元素 [a,b,c] ,s1=fa+fb+fc  增量为x，

-->[a+x,b+x,c+x] ,  s2=f(a+x)+f(b+x),f(c+x);

由矩阵关于斐波那契的递推式有:(fn,f_n-1)=(f_n-1,f_n-2)*(_{1 0}^{1 1})    //请自行脑补对齐= =

-->(f(a+x),f(a+x-1))=(f(a),f(a-1))*(_{1 0}^{1 1})^x

(f(b+x),f(b+x-1))=(f(b),f(b-1))*(_{1 0}^{1 1})^x

(f(c+x),f(c+x-1))=(f(c),f(c-1))*(_{1 0}^{1 1})^x

将三个等式相加得到 (s2,f(a+x-1)+f(b+x-1)+f(c+x-1))=(s1,f(a-1)+f(b-1)+f(c-1))*(_{1 0}^{1 1})^x  //由于矩阵运算满足单向的分配律

我们已经找到了s2和s1的关系了，在已知x和行矩阵的第二个元素的情况下，只需要对(_{1 0}^{1 1})进行一次矩阵幂便可得到s2.

由于计算时必须要知道行矩阵的第二项，我们不妨对于每个节点维护两个值，是Σ^r_i=l f(i) 和  Σ^r_i=l f(i-1)，利用laz标记的x便可快速的完成释放操作。

  1 #include<bits/stdc++.h>
  2 using namespace std;
  3 #define LL long long
  4 LL MOD=1e9+7;
  5 const LL MAXN=(100000<<2)+15;
  6 struct Matrix
  7 {
  8     int r=2,w=2;
  9     LL a[4][4];
 10     void init(){memset(a,0,sizeof(a));}
 11     Matrix operator*(const Matrix &tmp){
 12      Matrix ans;
 13      ans.r=r;ans.w=tmp.w;
 14      ans.init();
 15      for(int i=1;i<=r;++i)
 16         for(int k=1;k<=w;++k)
 17           for(int j=1;j<=tmp.w;++j)
 18     ans.a[i][j]=(ans.a[i][j]+a[i][k]*tmp.a[k][j])%MOD;
 19     return ans;
 20     }
 21     Matrix operator+(const Matrix &tmp){
 22      Matrix ans;
 23      ans.r=r;ans.w=w;
 24      ans.init();
 25      for(int i=1;i<=r;++i)
 26         for(int j=1;j<=w;++j)
 27     ans.a[i][j]=(ans.a[i][j]+a[i][j]+tmp.a[i][j])%MOD;
 28     return ans;
 29     }
 30 }unit,Am;
 31 void init()
 32 {
 33     unit.init();unit.r=unit.w=2;
 34     for(int i=0;i<=3;++i) unit.a[i][i]=1;
 35     Am.init(),Am.r=Am.w=2;
 36     Am.a[1][1]=Am.a[1][2]=Am.a[2][1]=1;
 37 }
 38 Matrix qpow(Matrix A,LL n)
 39 {
 40  Matrix ans=unit;
 41  ans.r=A.r;
 42  ans.w=A.w;
 43  while(n){
 44     if(n&1) ans=ans*A;
 45     A=A*A;
 46     n>>=1;
 47  }
 48  return ans;
 49 }
 50 struct node
 51 {
 52     LL x,y;
 53 };
 54 struct Segtree
 55 {
 56     #define M ((L+R)>>1)
 57     #define lc (id<<1)
 58     #define rc (id<<1|1)
 59     LL laz[MAXN];
 60     node sum[MAXN];
 61     void init()
 62     {
 63         memset(laz,0,sizeof(laz));
 64         memset(sum,0,sizeof(sum));
 65     }
 66     void build(int L,int R,int id)
 67     {
 68        if(L==R){
 69         LL x; scanf("%lld",&x);
 70         if(x==1){sum[id].x=1;}
 71         else if(x==2) {sum[id].x=sum[id].y=1;}
 72         else{
 73             Matrix t=qpow(Am,x-2);
 74             sum[id].x=(t.a[1][1]+t.a[1][2])%MOD;
 75             sum[id].y=(t.a[1][2]+t.a[2][2])%MOD;
 76            // cout<<sum[id].x<<" "<<sum[id].y<<endl;
 77         }
 78         return;
 79        }
 80        build(L,M,lc);
 81        build(M+1,R,rc);
 82        pushup(L,R,id);
 83     }
 84     void pushdown(int L,int R,int id)
 85     {
 86      if(!laz[id]) return;
 87      laz[lc]+=laz[id];
 88      laz[rc]+=laz[id];
 89      Matrix t=qpow(Am,laz[id]);
 90      LL x1=sum[lc].x,y1=sum[lc].y;
 91      sum[lc].x=(x1*t.a[1][1]%MOD+y1*t.a[2][1]%MOD)%MOD;
 92      sum[lc].y=(x1*t.a[1][2]%MOD+y1*t.a[2][2]%MOD)%MOD;
 93       x1=sum[rc].x,y1=sum[rc].y;
 94      sum[rc].x=(x1*t.a[1][1]%MOD+y1*t.a[2][1]%MOD)%MOD;
 95      sum[rc].y=(x1*t.a[1][2]%MOD+y1*t.a[2][2]%MOD)%MOD;
 96      laz[id]=0;
 97     }
 98     void pushup(int L,int R,int id)
 99     {
100      sum[id].x=(sum[lc].x+sum[rc].x)%MOD;
101      sum[id].y=(sum[lc].y+sum[rc].y)%MOD;
102     }
103     void update(int L,int R,int id,int l,int r,LL v)
104     {
105      if(L>=l&&R<=r){
106         laz[id]+=v;
107         Matrix t=qpow(Am,v);
108         LL x1=sum[id].x,y1=sum[id].y;
109         sum[id].x=(x1*t.a[1][1]%MOD+y1*t.a[2][1]%MOD)%MOD;
110         sum[id].y=(x1*t.a[1][2]%MOD+y1*t.a[2][2]%MOD)%MOD;
111         return;
112      }
113      pushdown(L,R,id);
114      if(l<=M) update(L,M,lc,l,r,v);
115      if(r>M) update(M+1,R,rc,l,r,v);
116      pushup(L,R,id);
117     }
118     LL ask(int L,int R,int id,int l,int r)
119     {
120      if(L>=l&&R<=r) return sum[id].x%MOD;
121      pushdown(L,R,id);
122      LL s=0;
123      if(l<=M) s=(s+ask(L,M,lc,l,r))%MOD;
124      if(r>M) s=(s+ask(M+1,R,rc,l,r))%MOD;
125      pushup(L,R,id);
126      return s;
127     }
128 }seg;
129 int main()
130 {
131     init();
132     LL n,m,x,a;
133     int i,opt,j,l,r,k;
134     while(cin>>n>>m){
135         seg.init();
136         seg.build(1,n,1);
137         for(i=1;i<=m;++i){
138             cin>>opt;
139             LL v;
140             if(opt==1){
141                 cin>>l>>r>>v;
142                 seg.update(1,n,1,l,r,v);
143             }
144             else{
145                 cin>>l>>r;
146                 cout<<seg.ask(1,n,1,l,r)<<endl;
147             }
148         }
149     }
150     return 0;
151 }
152  

1284: SP教数学

时间限制: 2 秒  内存限制: 128 MB
提交: 24  解决: 4

题目描述

输入

在

---恢复内容结束---