题解：CF226C Anniversary

题目传送门

说明

设 $F$ 为斐波那契数列。

性质

首先，你需要知道一个性质：
$$\gcd (F_n,F_{n+1})=1$$
证明可以用辗转相减法，由于本人不会 Latex，那么直接用文字证明：

• $F_{n+1}-Fn=Fn-1$，所以：
• $\gcd (F_n,F_{n+1})=\gcd (F_{n-1},F_n)$。
• 那么，当 $n$ 递归到 2 的时候就得 $\gcd (1,1)=1$。
• 固得证。

之后，我们要证明一个东西。

设 $n<m$，$x=F_n$，$y=F_{n+1}$，那么：

• $F_{n+1}=x+y$，$F_{n+2}=x+2y$，$F_{n+3}=2x+3y$……$F_m=F_{m-n-1}\times x+F_{m-n}\times b$。
• 那么，$F_m=F_{m-n-1}\times F_n+F_{m-n}\times F_{n+1}$。
• 那么，$\gcd (F_n,F_m)=\gcd (F_n,F_{m-n}\times F_{n+1})$。
• 由互质性质可知，$\gcd (F_n,F_m)=\gcd (F_n,F_{m-n})$。
• 递归下去，可以知道答案就是 $F_{\gcd (n,m)}$。

分析

我们需要求出一个最大的 $res$，使得 ${\sum }_{i=l}^r[res|F_i]\ge k$。由于上述的性质，所以可以变为：${\sum }_{i=1}^r[res|i]\ge k$。注意！这个 $res$ 不可以二分。证明：

• 设计数函数 $cnt(res)=\lfloor \frac{r}{res}\rfloor -\lfloor \frac{l-1}{res}\rfloor$。那么二分所需要的单调性是 $cnt(res)\ge cnt(res+1)$。
• 但是，$l=4,r=8$ 的情况可以知道 $cnt(3)<cnt(4)$，原因是向下取整会有精度流失，不能够强行证明。而且如果两端都是 $res+1$ 的倍数，那么 $cnt(res+1)\ge 2$，但是只能保证 $cnt(res)\ge 1$。

那么怎么样才能够求出这一个 $res$ 呢？其实可以用根号枚举。我们需要的是 $cnt(res)$ 不同的结果，这样的结果最多有 $\lfloor \sqrt{r}\rfloor$ 个，有时候会多一些。那么，我们求：
$$res=\underset{i=1}{\overset{\lfloor \sqrt{r}\rfloor }{\max }}(\max (cnt(i),cnt(\lfloor \frac{n}{i}\rfloor )))$$
即可。

加速

可以发现这一个值可能会很大，那么我们需要很高效的算法，求出 $F_{res}$。这个算法就是矩阵快速幂优化 dp。

可以发现斐波那契数列本质上是一个 dp，即：
$$d{p}_i=\{\begin{array}{ll}1& i\le 2\\ d{p}_{i-2}+d{p}_{i-1}& i>3\end{array}$$
矩阵快速幂的本质是利用矩阵乘法的乘法结合律来用 $\log n$ 的时间复杂度来求出分矩阵，然后乘上原矩阵。

矩阵乘法的规则如下：

• 定义 $x$ 为长度 $n,m$ 的矩阵，$y$ 为长度 $m,k$ 的矩阵，那么它们相乘的结果为长度为 $n,k$ 的矩阵 $res$。
• $re{s}_{i,j}={\sum }_{k=1}^m{x}_{i,k}\times y_{k,j}$。

可以看出这一个式子似乎能够优化所有转移方程只带加法和乘法的 dp。但是，要证明两点：

• 矩阵乘法不满足乘法交换律，因为交换后 $n$ 可能不等于 $k$。
• 矩阵乘法满足结合律。因为 { n , m } × { m , k } × { k , l } \{n,m\}\times\{m,k\}\times\{k,l\} {n,m}×{m,k}×{k,l}，不管先算哪一个，要么得出 { n , k } × { k , l } = { n , l } \{n,k\}\times\{k,l\}=\{n,l\} {n,k}×{k,l}={n,l}，要么得出 { n , m } × { m , l } = { n , l } \{n,m\}\times\{m,l\}=\{n,l\} {n,m}×{m,l}={n,l}。

所以，矩阵乘法满足结合律。那么，矩阵乘法的递推式要视情况而定。比如这一道题目需要：

$$\left(\begin{array}{cc}0& 1\end{array}\right)\times {\left(\begin{array}{cc}1& 1\\ 1& 0\end{array}\right)}^{res}$$
因为可以知道 $res$ 等于 1 的时候，乘完之后会得：
$$\left(\begin{array}{cc}1& 1\end{array}\right)$$
这时候可能不太明显，再乘 $n$ 次：
$$\left(\begin{array}{cc}1& 2\end{array}\right)$$
$$\left(\begin{array}{cc}2& 3\end{array}\right)$$
$$\left(\begin{array}{cc}3& 5\end{array}\right)$$
$$\left(\begin{array}{cc}5& 8\end{array}\right)$$
这是斐波那契数列。为什么是斐波那契数列呢？因为可以变成：
$$\left(\begin{array}{cc}d{p}_{i-1}& d{p}_i\end{array}\right)\times \left(\begin{array}{cc}1& 1\\ 1& 0\end{array}\right)=\left(\begin{array}{cc}d{p}_i& d{p}_{i-1}+d{p}_i\end{array}\right)$$
这就是矩阵快速幂优化 dp 的过程。在一些线性 dp 会超时的时候就可以用矩阵快速幂优化。可以优化到 $\log n$ 级别的递推式。

所以，可以构造如上的矩阵然后求矩阵快速幂。

代码

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
struct matrix{
	int n,m;
	ll a[3][3];
	matrix(){
		memset(a,0,sizeof(a));
	}
	inline ll *operator[](int index){
		return a[index];
	} 
};
int mod;
ll l,r,k;
inline matrix operator*(matrix x,matrix y){
	matrix res;
	res.n=x.n;
	res.m=y.m;
	for(int i=1;i<=x.n;++i){
		for(int j=1;j<=x.m;++j){
			for(int k=1;k<=y.m;++k){
				(res[i][k]+=x[i][j]*y[j][k])%=mod;
			}
		}
	}
	return res;
}
inline matrix power(matrix x,ll y){
	matrix res=x;
	while(y){
		if(y&1){
			res=res*x;
		}
		y>>=1;
		x=x*x;
	}
	return res;
}
int main(){
	scanf("%d %lld %lld %lld",&mod,&l,&r,&k);
	ll res=0;
	for(ll i=1;i*i<=r;++i){
		if(r/i-(l-1)/i>=k){
			res=max(res,i);
		}
		if(r/(r/i)-(l-1)/(r/i)>=k){
			res=max(res,r/i);
		}
	}
	matrix x,y;
	x.n=x.m=y.n=2;
	y.m=1;
	x[1][1]=x[1][2]=x[2][1]=y[2][1]=1;
	printf("%lld",matrix(power(x,res-1)*y)[1][1]);
	return 0;
}