概述

对于类似于 $(a+\sqrt{b})^n$ 的式子
我们可以共轭构造如下递推式

$设0 < a - \sqrt{b} < 1,令S_n = ((a+\sqrt{b})^n + (a-\sqrt{b})^n)$
$\begin{align}
&S_n((a+\sqrt{b}) + (a-\sqrt{b})) \\
&= ((a+\sqrt{b})^n + (a-\sqrt{b})^n)((a+\sqrt{b}) + (a-\sqrt{b})) \\
&= (a+\sqrt{b})^{n+1} + (a-\sqrt{b})^{n+1} + (a+\sqrt{b})^n(a-\sqrt{b}) + (a-\sqrt{b})^n(a+\sqrt{b}) \\
&= S_{n+1} + (a^2 - b)((a+\sqrt{b})^{n-1} + (a-\sqrt{b})^{n-1}) \\
&= S_{n+1} + (a^2 - b)S_{n-1} \\
&\Rightarrow S_{n+1} = 2aS_n - (a^2 - b)S_{n-1}
\end{align}
$

有了递推式,便可以用矩阵快速幂求解

题目

hdu4565

求$S_n = \lceil (a+\sqrt{b})^n \rceil \mod m$
题目中有一个很重要的条件$(a-1)^2 < b < a^2$,这样可以直接套用上面的递推式然后用矩阵快速幂求解

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
#include<cstdio>
#include<cstring>

const int maxn = 3;
typedef long long LL;
typedef LL Matrix[maxn][maxn];

const int sz = 2; //矩阵大小
LL mod; //模

void matrix_mul(Matrix A,Matrix B,Matrix res)
{
Matrix C;
memset(C,0,sizeof(C));
for(int i = 0;i < sz;++i)
for(int j = 0;j < sz;++j)
for(int k = 0;k < sz;++k)
C[i][j] = ((C[i][j] + A[i][k] * B[k][j])%mod + mod) % mod;
memcpy(res,C,sizeof(C));
}

void matrix_pow(Matrix A,LL n,Matrix res)
{
Matrix a,r;
memcpy(a,A,sizeof(a));
memset(r,0,sizeof r);

for(int i=0;i<sz;++i) r[i][i] = 1;
while(n)
{
if(n&1) matrix_mul(r,a,r);
n >>= 1;
matrix_mul(a,a,a);
}
memcpy(res,r,sizeof(r));
}

int main()
{
//freopen("test.txt","r",stdin);
LL a,b,n;
while(scanf("%lld %lld %lld %lld",&a,&b,&n,&mod)==4)
{
Matrix A,T,res;
memset(A,0,sizeof A);
memset(T,0,sizeof T);

T[0][0] = 2*(a*a +b)%mod ; // S2
T[1][0] = 2*a % mod; // S1

A[0][0] = 2*a % mod;
A[0][1] = ((b - a*a)%mod +mod) % mod;
A[1][0] = 1;

if(n==1)
{
printf("%lld\n",T[1][0]);
continue;
}
else if(n==2)
{
printf("%lld\n",T[0][0]);
continue;
}

matrix_pow(A,n-2,res);
matrix_mul(res,T,res);

printf("%lld\n",((res[0][0]%mod)+mod) % mod);
}
return 0;
}