- 【题目链接】
-
【前置知识】
-
矩阵快速幂,矩阵加速
-
扩展欧拉定理
-
-
【解题思路】
看到这个数据范围,我们应该能想到,这是一道结论题。
再看一下SPOJ几页的 0.00 s 0.00\operatorname s 0.00s,我们应该会觉得这题正解 Θ ( 1 ) \Theta(1) Θ(1)。
那么,我们推一下结论吧。
首先,我们设 f ( n ) = 1 f(n)=1 f(n)=1,那么计算 F n , f F_{n,f} Fn,f 中调用 f f f 的次数即为 F n , f F_{n,f} Fn,f 的返回值。(这里 F F F 请用题目翻译中的定义。)
然后,我们又可以知道, F a , F b , f = F a , f F b , f F_{a,F_{b,f}}=F_{a,f}F_{b,f} Fa,Fb,f=Fa,fFb,f。因为,显然, F a , F b , f F_{a,F_{b,f}} Fa,Fb,f 会把 F b , f F_{b,f} Fb,f 调用 F a , f F_{a,f} Fa,f 次。
看起来也没有什么思路,我最开始用的是暴力展开。然后,有一个问题,就是不可避免的要计算任意 n n n 对应的 F n , f F_{n,f} Fn,f 的值。于是又绕回了原来的问题。不过,我们可以试着,如果把所有的 F n , f F_{n,f} Fn,f 的值列出来,会怎样呢?
接着,我们就来从 F 1 , f F_1,f F1,f 开始写出它对应的值:
F 1 , f = a , F 2 , f = b , F 3 , f = F 2 , F 1 , f = a b , F 4 , f = F 3 , F 2 , f = a b 2 , F 5 , f = F 4 , F 3 , f = a 2 b 3 ⋯ F_{1,f}=a,F_{2,f}=b,F_{3,f}=F_{2,F_{1,f}}=ab,F_{4,f}=F_{3,F_{2,f}}=ab^2,F_{5,f}=F_{4,F_{3,f}}=a^2b^3\cdots F1,f=a,F2,f=b,F3,f=F2,F1,f=ab,F4,f=F3,F2,f=ab2,F5,f=F4,F3,f=a2b3⋯
发现了吗? a a a 与 b b b 的指数是斐波那契数列!
其实这个性质的证明也很显然,因为把前两项指数相加就是这一项的指数,而这正是斐波那契数列的定义。
那么,我们对斐波那契数列的边界条件做点微调:设数列 g g g,它满足以下条件:
g i = { 1 i = 1 0 i = 2 g i − 1 + g i − 2 i > 2 g_i= \begin{cases} 1&i=1\\ 0&i=2\\ g_{i-1}+g_{i-2}&i>2 \end{cases} gi=⎩⎪⎨⎪⎧10gi−1+gi−2i=1i=2i>2
则答案即为 a g n b g n + 1 m o d p a^{g_n}b^{g_{n+1}}\bmod p agnbgn+1modp。
而我们发现, g n g_n gn 与 g n + 1 g_{n+1} gn+1 可以通过矩阵加速 Θ ( log 2 n ) \Theta(\log_2 n) Θ(log2n) 求出!
那么,接下来就是一个快速幂,没了。复杂度 Θ ( log 2 n ) \Theta(\log_2n) Θ(log2n)。
没了?真没了?
好好再想一想,我们该对指数怎么处理?如果像平常那样把指数对 p p p 取模,你会发现样例都过不了。我们总不能用高精算一遍再取模吧?
你可能想,对 p − 1 p-1 p−1 取模不就行了?然后,你会发现,题目中没有说过 p p p 是质数,所以费马小定理不能用。 而且,最阴间的是,题目同样没有对 a , b a,b a,b 作额外限定。
那么,还有一种方法:扩展欧拉定理。它不需要 a , b a,b a,b 与 p p p 互质。
到这里,这个问题就被真正的解决了。
也许你可能会说:开始不是说要 Θ ( 1 ) \Theta(1) Θ(1) 吗?
是的。但是,这个 O ( p log 2 n ) \operatorname O(\sqrt p\log_2 n) O(plog2n) 算法确实是本题的正解。之所以我会说是 Θ ( 1 ) \Theta(1) Θ(1) 的结论题,是因为这个猜想会引导我们往这个方向去想,从而导出真正的正解。
- 【实现细节】
问题虽然被解决了,代码还是要注意的。
扩展欧拉定理怎么实现,就看你自己了。我的做法是在矩阵乘法函数里判定是否大于 φ ( p ) \varphi(p) φ(p)。
还有一种极为阴间的情况: p = 1 p=1 p=1!注意这一点,特判好了,就能过了。
对了,我曾经特判输出 0 0 0 的时候忘记输出测试信息了,调了好久……
- 【代码实现】
#include <iostream>
#include <vector>
using namespace std;
typedef unsigned long long lrg_l;
typedef vector<vector<lrg_l>> matrix;
lrg_l moder,prime;
matrix greater_varphi;
long long fast_pow(long long base,long long exp)
{
long long result;
for(result=1;exp;exp&1?result=result*base%prime:true,base=base*base%prime,exp>>=1);
return result;
}
matrix operator*(matrix& src,matrix& param)
{
const int len=src.size();
matrix result(len,vector<lrg_l>(len,0));
for(int k=0;k<len;k++)
for(int i=0;i<len;i++)
for(int j=0;j<len;j++)
{
result[i][j]=result[i][j]+src[i][k]*param[k][j];
if(result[i][j]>moder)
result[i][j]%=moder,greater_varphi[i][j]=true;
}
return result;
}
int main(int argc,char *argv[],char *envp[])
{
ios::sync_with_stdio(false);
cin.tie(nullptr);
long long testcnt,power,base_A,base_B,tmp;
cin>>testcnt;
matrix result(2,vector<lrg_l>(2)),base(2,vector<lrg_l>(2)),
initstate={{1,1},{1,0}},resstate={{0,1},{0,0}},greater_varphi_init={{false,false},{0,0}};
for(int t=1;t<=testcnt;t++)
{
cin>>base_A>>base_B>>prime>>power;
power--;
if(prime==1)
{
cout<<"Case #"<<t<<": 0\n";
continue;
}
base=initstate;
result=resstate;
greater_varphi=greater_varphi_init;
moder=prime;
tmp=prime;
for(int i=2;i*i<=prime;i++)
if(prime%i==0)
{
moder=moder/i*(i-1);
while(prime%i==0)
prime/=i;
}
if(prime>1)
moder=moder/prime*(prime-1);
prime=tmp;
while(power)
{
if(power&1)
result=result*base;
base=base*base;
power>>=1;
}
cout<<"Case #"<<t<<": "<<fast_pow(base_A,result[0][1]+greater_varphi[0][1]*moder)*fast_pow(base_B,result[0][0]+greater_varphi[0][0]*moder)%prime<<'\n';
}
return 0;
}