Luogu 块速递推 | 光速幂模板

Luogu P5110 块速递推

---

• 题意
  给定一个数列 (a) 满足递推式：
  [a_n=233a_{n-1}+666a_{n-2},a_0=0,a_1=1 ]
  求这个数列第 (n) 项 (mod) (10^9+7) 的值，有 (T) 组询问。
  (n) 在 unsigned long long 的范围内； (1 leq T leq 5*10^7)。

---

• 题解
  我们先列出这个递推式的特征方程：
  [233x+666=x^2 ]
  [x^2-233x-666=0 ]
  [x_1=frac{233+sqrt{56953}}{2};x_2=frac{233-sqrt{56953}}{2} ]
  [a_n=k_1x_1^n+k_2x_2^n ]
  [egin{cases} k_1+k_2=0 \ k_1x_1+k_2x_2=1 end{cases}]
  [egin{cases} k_1=frac{1}{sqrt{56953}} \ k_1=frac{1}{-sqrt{56953}} end{cases}]
  [ herefore a_n=frac{1}{sqrt{56953}}((frac{233+sqrt{56953}}{2})^n-(frac{233-sqrt{56953}}{2})^n) ]
  [ecause 188305837=sqrt{56953}(mod　10^9+7) ]
  [ecause a_n=233230706*(94153035^n-905847205^n) ]
  然后我们接下来要做的就是 (O(1)) 求这个式子（光速幂）
  • 我们可以设 (x^n=32769k+r)
    [f_1(n)=x^{32768n};f_2(n)=x^n ]
    [ herefore x^p=f_1(frac{p}{32768})f_2(p　mod　32768) ]

---

• 接下来就是一点也不激动人心的代码了：

#include<bits/stdc++.h>
/* Forever_chen */
#define RT register
#define int long long
#define P 1000000007
using namespace std;
template<class t> inline t read(t &x){
	char c=getchar();bool f=0;x=0;
	while(!isdigit(c)) f|=c=='-',c=getchar();
	while(isdigit(c))x=(x<<1)+(x<<3)+(c^48),c=getchar();
	if(f)x=-x;return x;
}
template<class t>inline void write(t x){
	if(x<0)putchar('-'),write(-x);
	else{if(x>9)write(x/10);putchar('0'+x%10);}
}
template<class t>inline void writeln(t x){
	write(x);putchar('
');
	return;
}
template<class t>inline void write_blank(t x){
	write(x);putchar(' ');
	return;
}
namespace Mker{
    unsigned long long SA, SB, SC;
    void init() { 
		scanf("%llu%llu%llu", &SA, &SB, &SC); 
	}
    inline unsigned long long rand(){
        SA ^= SA << 32, SA ^= SA >> 13, SA ^= SA << 1;
        unsigned long long t = SA;
        SA = SB, SB = SC, SC ^= t ^ SA;
		return SC;
    }
}
in
int x1,x2,x3,x4;
int f1[65550],f2[65550],f3[65550],f4[65550];
inline int Pow_1(int x) { 
	return 1ll*f3[x>>16]*f1[x&65535]%P; 
}
inline int Pow_2(int x) { 
	return 1ll*f4[x>>16]*f2[x&65535]%P; 
}
int ttt,ans;
signed main(){
	//freopen(".in","r",stdin);
	//freopen(".out","w",stdout);
	f1[0]=1;
	f2[0]=1;
	f3[0]=1;
	f4[0]=1;
	x1=94153035;
	x2=905847205; 
	x3=64353223;
	x4=847809841;
    for(RT int i=1;i<65536;i++){
    	f1[i]=1ll*f1[i-1]*x1%P;
	} 
    for(RT int i=1;i<65536;i++){
    	f2[i]=1ll*f2[i-1]*x2%P;
	} 
    for(RT int i=1;i<65536;i++){
    	f3[i]=1ll*f3[i-1]*x3%P;
	}
    for(RT int i=1;i<36636;i++){
    	f4[i]=1ll*f4[i-1]*x4%P;
	} 
	read(ttt);
	Mker::init();
	while(ttt--){
		unsigned int lop;
		lop=Mker::rand()%(P-1);
		ans^=233230706*(Pow_1(lop)-Pow_2(lop)+P)%P; 
	}
	write(ans);
	//system("pause");
	return 0;
}

注：可能需要卡卡常