(f_n=f_{n-1}+f_{n-2},f_1=f_2=1),求(f_1+2f_2+3f_3+...+nf_nmod m,1≤n,m≤2^{31}-1)。
解
数列问题加比较大的数据范围，就很容易到与矩阵快速幂有关，于是尝试变换式子,注意任何小的看起来不起眼的式子变换都有不同的结果,注意递推转移常用的不是策略，而是问题的划分
法一：
设(t_n=f_1+2f_2+...+nf_n),有(t_n=t_{n-1}+nf_n),现在关键在于求(nf_n),根据什么变维护什么的理论，设(g_n=nf_n=n(f_{n-1}+f_{n-2})=)
((n-1)f_{n-1}+(n-2)f_{n-2}+f_{n-1}+2f_{n-2}=g_{n-1}+g_{n-2}+f_{n-1}+2f_{n-2}),于是我们要想得到t，三个递推都得维护，所以不难有状态矩阵
[egin{bmatrix}f_{n-2}&f_{n-1}&g_{n-2}&g_{n-1}&t_{n}end{bmatrix}
]
按照填矩阵转移方程套路，不难有转移矩阵
[egin{bmatrix}0&1&0&2&0\1&1&0&1&0\0&0&0&1&0\0&0&1&1&1\0&0&0&0&1end{bmatrix}
]
按照基本套路转移即可。
参考代码：
#include <iostream>
#include <cstdio>
#include <cstring>
#define il inline
#define ri register
#define ll long long
using namespace std;
ll yyb;
struct matrix{
ll jz[5][5];
il void clear(){
memset(jz,0,sizeof(jz));
}
il void unit(){
clear();
for(ri int i(0);i<5;++i)jz[i][i]|=true;
}
il matrix operator*(matrix x){
matrix y;y.clear();
ri int i,j,k;
for(i=0;i<5;++i)
for(j=0;j<5;++j)
for(k=0;k<5;++k)
(y.jz[i][j]+=jz[i][k]*x.jz[k][j])%=yyb;
return y;
}template<class free>
il matrix operator^(free y){
matrix ans,x(*this);ans.unit();
while(y){
if(y&1)ans=ans*x;
x=x*x,y>>=1;
}return ans;
}
}s,t;
int main(){
s.jz[0][0]=0,s.jz[0][1]=1,s.jz[0][2]=0,s.jz[0][3]=1,s.jz[0][4]=0;
t.jz[0][0]=0,t.jz[0][1]=1,t.jz[0][2]=0,t.jz[0][3]=2,t.jz[0][4]=0;
t.jz[1][0]=1,t.jz[1][1]=1,t.jz[1][2]=0,t.jz[1][3]=1,t.jz[1][4]=0;
t.jz[2][0]=0,t.jz[2][1]=0,t.jz[2][2]=0,t.jz[2][3]=1,t.jz[2][4]=0;
t.jz[3][0]=0,t.jz[3][1]=0,t.jz[3][2]=1,t.jz[3][3]=1,t.jz[3][4]=1;
t.jz[4][0]=0,t.jz[4][1]=0,t.jz[4][2]=0,t.jz[4][3]=0,t.jz[4][4]=1;
ll n;scanf("%lld%lld",&n,&yyb),s=s*(t^n),printf("%lld",s.jz[0][4]);
return 0;
}
法二：
设(s_n=sum_{i=1}^nf_i)
[f_1+2f_2+3f_3+...+nf_n=s_n+s_n-s_1+s_n-s_2+...s_n-s_{n-1}
]
[=ns_n-(s_1+s_2+...+s_{n-1})
]
前面一截已经可以很好算了，于是考虑变换后面一截，单独拿出来考虑，设
[g(n)=s_1+s_2+...s_n=g_{n-1}+s_n
]
[s_n=s_{n-1}+f_n
]
[f_n=f_{n-1}+f_{n-2}
]
而原式为
[ans=ns_n-g_{n-1}
]
所以只要想办法求出这两个东西即可，于是考虑同时转移，所以设状态矩阵
[egin{bmatrix}f_n&f_{n+1}&s_n&g_{n-1}end{bmatrix}
]
根据填转移矩阵套路，不难有转移矩阵
[egin{bmatrix}0&1&0&0\1&1&1&0\0&0&1&1\0&0&0&1end{bmatrix}
]
以此转移即可,其实实质就是前缀和套前缀和再套前缀和，解释一下，s为f的前缀和，g是s的前缀和。
参考代码：
#include <iostream>
#include <cstdio>
#include <cstring>
#define il inline
#define ri register
#define ll long long
using namespace std;
ll yyb;
struct matrix{
ll jz[4][4];
il void clear(){
memset(jz,0,sizeof(jz));
}
il void unit(){
clear();
for(ri int i(0);i<4;++i)jz[i][i]=1;
}
il matrix operator*(matrix x){
matrix y;y.clear();ri int i,j,k;
for(i=0;i<4;++i)
for(j=0;j<4;++j)
for(k=0;k<4;++k)
(y.jz[i][j]+=jz[i][k]*x.jz[k][j]%yyb)%=yyb;
return y;
}template<class free>
il matrix operator^(free y){
matrix ans,x(*this);ans.unit();
while(y){
if(y&1)ans=ans*x;
x=x*x,y>>=1;
}return ans;
}
}state,tran;
int main(){
ll n,ans;
scanf("%lld%lld",&n,&yyb);
state.jz[0][0]=1,state.jz[0][1]=1,state.jz[0][2]=1,state.jz[0][3]=0;
tran.jz[0][0]=0,tran.jz[0][1]=1,tran.jz[0][2]=0,tran.jz[0][3]=0;
tran.jz[1][0]=1,tran.jz[1][1]=1,tran.jz[1][2]=1,tran.jz[1][3]=0;
tran.jz[2][0]=0,tran.jz[2][1]=0,tran.jz[2][2]=1,tran.jz[2][3]=1;
tran.jz[3][0]=0,tran.jz[3][1]=0,tran.jz[3][2]=0,tran.jz[3][3]=1;
state=state*(tran^n-1),ans=((state.jz[0][2]*n%yyb-state.jz[0][3])%yyb+yyb)%yyb;
printf("%lld",ans);
return 0;
}
#include <iostream>
#include <cstdio>
#include <cstring>
#define il inline
#define ri register
#define ll long long
using namespace std;
ll yyb;
struct matrix{
ll jz[4][4];
il void clear(){
memset(jz,0,sizeof(jz));
}
il void unit(){
clear();
for(ri int i(0);i<4;++i)jz[i][i]=1;
}
il matrix operator*(matrix x){
matrix y;y.clear();ri int i,j,k;
for(i=0;i<4;++i)
for(j=0;j<4;++j)
for(k=0;k<4;++k)
(y.jz[i][j]+=jz[i][k]*x.jz[k][j]%yyb)%=yyb;
return y;
}template<class free>
il matrix operator^(free y){
matrix ans,x(*this);ans.unit();
while(y){
if(y&1)ans=ans*x;
x=x*x,y>>=1;
}return ans;
}
}state,tran;
int main(){
ll n,ans;
scanf("%lld%lld",&n,&yyb);
state.jz[0][0]=1,state.jz[0][1]=1,state.jz[0][2]=1,state.jz[0][3]=0;
tran.jz[0][0]=0,tran.jz[0][1]=1,tran.jz[0][2]=0,tran.jz[0][3]=0;
tran.jz[1][0]=1,tran.jz[1][1]=1,tran.jz[1][2]=1,tran.jz[1][3]=0;
tran.jz[2][0]=0,tran.jz[2][1]=0,tran.jz[2][2]=1,tran.jz[2][3]=1;
tran.jz[3][0]=0,tran.jz[3][1]=0,tran.jz[3][2]=0,tran.jz[3][3]=1;
state=state*(tran^n-1),ans=((state.jz[0][2]*n%yyb-state.jz[0][3])%yyb+yyb)%yyb;
printf("%lld",ans);
return 0;
}
法三：
实际上我们可能有些东西不需要转移，递推里面很多式子都会有通向公式，而我们有结论(s_n=sum_{i=1}^nf_i=f_{n+2}-f_2)，接着想办法优化，接着法二
[ans=ns_n-(s_1+s_2+...+s_{n-1})=
]
[ns_n-(f_3-f_2+...+f_{n+1}-f_2)=ns_n-(s_{n+1}-n-1)
]
[=n(f_{n+2}-f_2)-(f_{n+3}-n-2)=
]
[=nf_{n+2}-f_{n+3}+n+2
]
于是我们只要递推处f就可以算出ans了。
参考代码：
#include <iostream>
#include <cstdio>
#include <cstring>
#define il inline
#define ri register
#define ll long long
using namespace std;
ll yyb;
struct matrix{
ll jz[2][2];
il void clear(){
memset(jz,0,sizeof(jz));
}
il void unit(){
clear();
for(ri int i(0);i<2;++i)jz[i][i]|=true;
}
il matrix operator*(matrix x){
matrix y;y.clear();
ri int i,j,k;
for(i=0;i<2;++i)
for(j=0;j<2;++j)
for(k=0;k<2;++k)
(y.jz[i][j]+=jz[i][k]*x.jz[k][j]%yyb)%=yyb;
return y;
}template<class free>
il matrix operator^(free y){
matrix ans,x(*this);ans.unit();
while(y){
if(y&1)ans=ans*x;
x=x*x,y>>=1;
}return ans;
}
}state,tran;
int main(){
ll n,ans;
scanf("%lld%lld",&n,&yyb);
state.jz[0][0]=0,state.jz[0][1]=1;
tran.jz[0][0]=0,tran.jz[0][1]=1;
tran.jz[1][0]=1,tran.jz[1][1]=1;
state=state*(tran^n+2),ans=state.jz[0][0]*n-state.jz[0][1]+2;
printf("%lld",(ans%yyb+yyb)%yyb);
return 0;
}
小结
不难看出随着优化程度的提高，问题解决办法也就越来越间接，矩阵维数也就越来越少，而这个优化的关键在于通项。