【2016 ACM-ICPC 曼谷区域赛 Gym-101161 G】Binary Strings【矩阵快速幂】

题意：

$01$串，不包含两个连续的 $1$。要求给出满足下述三个条件的 $01$ 串的个数。(答案模 $1e9+7$)

• $01$串长度在 $[L,R]$ 之间 $(1\le L\le R\le 10^{18})$
• $01$串长度是 $K$ 的倍数 $(1\le K\le 10^9)$
• 不包含连续两个$1$的$01$串

---

思路：

$dp[i][0/1]$ 表示长度为$i$，末尾为$0/1$的字符串个数。
$dp[i+1][0]=dp[i][0]+dp[i][1],dp[i+1][1]=dp[i][0]$。

通过打表，可以发现 $dp[i]$ 是一个斐波那契数列，$dp[1]=2,dp[2]=3,dp[3]=5,dp[n]=dp[n-1]+dp[n-2]$。(也可以直接推导得到该性质)

因此先通过斐波那契数列的常规求解方式求出 $dp[k]$ 与 $dp[k-1]$。再求出 $sum[x]$ ，$sum[x]$ 表示长度小于等于$x$，且长度是 $K$ 倍数的 $01$ 串总数。

假设 $a$ 是最大的小于 $L$ 的 $K$ 的倍数，$b$ 是最大的小于 $R$ 且大于 $L$ 的 $K$ 的倍数，则 $ans=sum[b]-sum[a]$。

$$\left(\begin{array}{ccc}dp[i-1]& dp[i]& sum[i]\end{array}\right)\ast \left(\begin{array}{ccc}0& 1& 1\\ 1& 1& 1\\ 0& 0& 1\end{array}\right)=\left(\begin{array}{ccc}dp[i]& dp[i+1]& sum[i+1]\end{array}\right)$$

但是 $sum[x]$ 只表示长度为 $K$ 倍数的 $01$ 串的总数，因此仅当 $dp[i]$ 的 $i$ 为 $K$ 的倍数时，$dp[i]$ 才会被累加到答案 $sum[i]$ 中。

$$\left(\begin{array}{ccc}dp[k-1]& dp[k]& sum[k]\end{array}\right)\ast {\left(\begin{array}{ccc}0& 1& 0\\ 1& 1& 0\\ 0& 0& 1\end{array}\right)}^{k-1}\ast \left(\begin{array}{ccc}0& 1& 1\\ 1& 1& 1\\ 0& 0& 1\end{array}\right)=\left(\begin{array}{ccc}dp[2\ast k-1]& dp[2\ast k]& sum[2\ast k]\end{array}\right)$$

此处 $sum[2\ast k]=sum[k]+dp[2\ast k]$，因此先矩阵快速幂求出 $k-1$ 阶的矩阵，乘上后面的 $1$ 阶矩阵，得到矩阵 $base$，再对 $base$ 进行快速幂，求出 $sum[a]$ 与 $sum[b]$，即可得到答案。

---

总结：

此题的主要难点应该在于打表，如果不能发现 $dp$ 数组是斐波那契数组的话，还是非常难以进行后续矩阵快速幂的推导的。

注意矩阵快速幂的矩阵乘积形式一定要写正确，比赛时将左边的横矩阵写成了列矩阵，导致一直 $Wa$，很是痛心！

---

代码：

代码改来改去的，非常丑陋，不具有借鉴价值…

#include <cstdio>
#include <iostream>
#include <cstring>
#define rep(i,a,b) for(int i = a; i <= b; i++)
#define LOG1(x1,x2) cout << x1 << ": " << x2 << endl;
#define LOG2(x1,x2,y1,y2) cout << x1 << ": " << x2 << " , " << y1 << ": " << y2 << endl;
#define LOG3(x1,x2,y1,y2,z1,z2) cout << x1 << ": " << x2 << " , " << y1 << ": " << y2 << " , " << z1 << ": " << z2 << endl;
typedef long long ll;
using namespace std; 
const int mod = 1e9+7;
typedef long long ll;
 
ll L,R,K;
//i、k、j只能作为临时变量 
struct Matrix{
	ll a[10][10];
	Matrix(){memset(a,0,sizeof(a));}
	
	Matrix operator * (const Matrix y)
	{
		Matrix ans;
		for(int i = 1;i <= 3;i++)   //行 
			for(int j = 1;j <= 3;j++)  //列 
				for(int k = 1;k <= 3;k++)
					ans.a[i][j] = (ans.a[i][j]+a[i][k]*y.a[k][j]%mod)%mod;			
		return ans;
	}
};
 
Matrix q_pow(Matrix x,ll k)
{
	Matrix ans = x; k--;
	while(k > 0)
	{
		if(k&1)
			ans = ans*x;
		x = x*x;
		k = k>>1;
	}
	return ans;
}
 
ll solve(ll a1,ll a2,ll sum1, Matrix x,ll n){ //求sum[k], n >= 2
	Matrix tp = q_pow(x,n-1ll);
	ll ans = ((tp.a[1][3]*a1%mod+tp.a[2][3]*a2%mod)%mod+tp.a[3][3]*sum1%mod)%mod;
	return ans;
}
 
int main()
{
	int _; scanf("%d",&_);
	rep(Ca,1,_)
	{
		scanf("%lld%lld%lld",&L,&R,&K);
		Matrix p;
		for(int i = 1;i <= 3;i++)
			for(int j = 1;j <= 3;j++)
				p.a[i][j] = 1;
		p.a[1][1] = p.a[3][1] = p.a[3][2] = 0;
		ll a1 = 2, a2 = 3, sum1 = 5;
		if(K >= 3){
			Matrix tmp = q_pow(p,K-2ll);
			a1 = ((tmp.a[1][1]*2ll%mod+tmp.a[2][1]*3ll%mod)%mod+tmp.a[3][1]*5ll%mod)%mod;
			a2 = ((tmp.a[1][2]*2ll%mod+tmp.a[2][2]*3ll%mod)%mod+tmp.a[3][2]*5ll%mod)%mod;
			sum1 = a2;
		}
		else if(K == 1){
			ll ans1 = 0, ans2 = 0;
			if(L == 1) ans1 = 0;
			else if(L == 2) ans1 = 2;
			else if(L == 3) ans1 = 5;
			else{
				Matrix tmp = q_pow(p,L-3);
				ans1 = ((tmp.a[1][3]*2ll%mod+tmp.a[2][3]*3ll%mod)%mod+tmp.a[3][3]*5ll%mod)%mod;
			}
 
			if(R == 1) ans2 = 2;
			else if(R == 2) ans2 = 5;
			else{
				Matrix tmp = q_pow(p,R-2);
				ans2 = ((tmp.a[1][3]*2ll%mod+tmp.a[2][3]*3ll%mod)%mod+tmp.a[3][3]*5ll%mod)%mod;
			}
			printf("Case %d: %lld\n",Ca,(ans2-ans1+mod)%mod);
			continue;
		}
		else sum1 = 3;
 
		Matrix tp1,tp2;
		for(int i = 1;i <= 3;i++)
			for(int j = 1;j <= 3;j++)
				tp1.a[i][j] = 1;
		tp1.a[1][1] = tp1.a[3][1] = tp1.a[3][2] = tp1.a[1][3] = tp1.a[2][3] = 0ll;
 
		for(int i = 1;i <= 3;i++)
			for(int j = 1;j <= 3;j++)
				tp2.a[i][j] = 1;
		tp2.a[1][1] = tp2.a[3][1] = tp2.a[3][2] = 0;
 
		Matrix base = q_pow(tp1,K-1ll)*tp2;	
		ll hp1,hp2;
		hp1 = (L/K)*K;
		if(hp1 == L) hp1 = (L/K-1ll)*K;
 
		hp2 = R/K;
		ll ans1 = 0ll,ans2 = 0ll;
		if(hp2 > 1ll){
			ans1 = solve(a1,a2,sum1,base,hp2);
		}
		else if(hp2 == 1) ans1 = sum1;
 
		if(hp1 > K) ans2 = solve(a1,a2,sum1,base,hp1/K);
		else if(hp1 == K) ans2 = sum1;
		// LOG2("ans1",ans1,"ans2",ans2);
		printf("Case %d: %lld\n",Ca,(ans1-ans2+mod)%mod);
	}
	return 0;
}