HDU 6395 sequence （矩阵快速幂+区间分块+细节）*

#include<iostream>
#include<algorithm>
#include<string>
#include<map>//int dx[4]={0,0,-1,1};int dy[4]={-1,1,0,0};
#include<set>//int gcd(int a,int b){return b?gcd(b,a%b):a;}
#include<vector>
#include<cmath>
#include<stack>
#include<string.h>
#include<stdlib.h>
#include<cstdio>
#define ll long long
#define MAX 1000000000
#define ms memset
#define maxn 100005
#define INF 1000000
using namespace std;#define debug(x)             for(int i=0;i<3;i++,puts(""))  for(int j=0;j<3;j++)   cout<<x.a[i][j]<<" ";/*
题目意思很简单明了。
比赛菜的没做出来，赛后看标程才明白自己这部分已经生疏了很多，
记得以前有做过利用p/i除数相同的性质来分块的题目，好像是紫书上的题目，
利用这个求余数和的。
这题也是利用这个性质，比如p=16,
那么：i按照p/(p/i)+1来递增，
得到：1，2，3，4，5，6，9，17.
滑动的区间中每个数对于p来说除数一样，比如【9，17）这部分对于16来说除数都是1，那么观察如果后面的p/n是个常数，就是赤裸裸的矩阵快速幂了，
如果加上p/n,就要利用分块的思想，思路大体是这个思路，但实现起来有超多细节考虑
（真的不能理解赛场上为什么别人做的那么快，，或许是真的码少了）n与p的大小关系要界定，
n>=p时候，不仅要进行n<p的情况，最后要对n-p+1进行矩阵快速幂。
在循环的时候，i以p/（p/i）+1的形式推进。*/ll a,b,c,d,p,n;const int mod=1e9+7;
struct node{ll  a[3][3];node(){memset(a,0,sizeof(a));}node operator*(const node &y)///矩阵快速幂模板{node ans;for(int i=0;i<3;i++)for(int j=0;j<3;j++){ans.a[i][j]=0;for(int k=0;k<3;k++)   ans.a[i][j]=(ans.a[i][j]+a[i][k]*y.a[k][j]%mod)%mod;}return ans;}
};node ans,unit;
node matrix(node shu,int k){node ret;for(ret=unit;k;k>>=1,shu=shu*shu)  if(k&1) ret=ret*shu;return ret;
}node mi;
void init(ll c,ll d,ll e)
{memset(mi.a,0,sizeof(mi.a));mi.a[0][1]=1;mi.a[1][0]=c;mi.a[1][1]=d;mi.a[1][2]=e;mi.a[2][2]=1;
}int main()
{for(int i=0;i<3;i++) unit.a[i][i]=1;///单位矩阵int t;scanf("%d",&t);while(t--){scanf("%lld%lld%lld%lld%lld%lld",&a,&b,&c,&d,&p,&n);if(n==1) {printf("%lld\n",a);continue;}else if(n<p){memset(ans.a,0,sizeof(ans.a));ans.a[0][0]=a,ans.a[1][0]=b,ans.a[2][0]=1;for(ll i=3;i<=n;i=p/(p/i)+1){init(c,d,p/i);if(n>p/(p/i))    mi=matrix(mi,p/(p/i)+1-i);else mi=matrix(mi,n-i+1);ans=mi*ans;}printf("%lld\n", ans.a[1][0]);}else if(n>=p){memset(ans.a,0,sizeof(ans.a));ans.a[0][0]=a,ans.a[1][0]=b,ans.a[2][0]=1;for(ll i=3;i<=p;i=p/(p/i)+1){init(c,d,p/i);if(n>p/(p/i))    mi=matrix(mi,(p/(p/i))+1-i);else mi=matrix(mi,n-i+1);ans=mi*ans;}init(c,d,0);mi.a[1][2]=0;mi=matrix(mi,n-max(p,2LL));ans=mi*ans;printf("%lld\n",ans.a[1][0]);}}return 0;
}