首先得出,
b+d√2
和
b−d√2
是一元二次方程
x2−bx+b2−d4=0
的两根。
把
x2−bx+b2−d4=0
移项得
x2=bx+d−b24
。
两边同乘以
xn−2
,可以得出,
xn=bxn−1+d−b24xn−2
。
设
f[i]=(b+d√2)i+(b−d√2)i
,
此时就容易得出
f[i]
是个整数,并且递推式为
f[i]=bf[i−1]+d−b24f[i−2]
,
f[0]=2,f[1]=b
。这时候就能通过矩阵乘法求得
f[n]
(注意,相乘会爆long long,因此要用快速乘)。
最后考虑怎样通过
f[n]
求得结果。由于题目限定
b2≤d<(b+1)2
,所以
n
是奇数时
注意特判
n=0
时结果为
1
<script type="math/tex" id="MathJax-Element-24">1</script>。
代码:
#include <cmath>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
using namespace std;
typedef long long ll;
ll b, d, n, tm; const ll ZZQ = 7528443412579576937ll;
ll add(ll a, ll b) {
return (1ull * a + 1ull * b) % ZZQ;
}
ll prod(ll a, ll b) {
ll res = 0;
while (b) {
if (b & 1) res = add(res, a);
a = add(a, a);
b >>= 1;
}
return res;
}
struct cyx {
int n, m; ll v[4][4];
cyx() {}
cyx(int _n, int _m) :
n(_n), m(_m) {memset(v, 0, sizeof(v));}
friend inline cyx operator * (cyx a, cyx b) {
int i, j, k; cyx res = cyx(a.n, b.m);
for (i = 1; i <= res.n; i++) for (j = 1; j <= res.m; j++)
for (k = 1; k <= a.m; k++)
res.v[i][j] = add(res.v[i][j], prod(a.v[i][k], b.v[k][j]));
return res;
}
friend inline cyx operator ^ (cyx a, ll b) {
int i; cyx res = cyx(a.n, a.m);
for (i = 1; i <= res.n; i++) res.v[i][i] = 1;
while (b) {
if (b & 1) res = res * a;
a = a * a;
b >>= 1;
}
return res;
}
} P, Q;
int main() {
cin >> b >> d >> n; P = cyx(2, 2); Q = cyx(2, 1);
if (!n) return printf("1\n"), 0;
tm = (d >> 2) - prod(b + 1 >> 1, b - 1 >> 1);
P.v[1][1] = b; P.v[1][2] = tm; P.v[2][1] = 1;
Q.v[1][1] = b; Q.v[2][1] = 2; P = (P ^ n - 1) * Q;
ll ans = P.v[1][1]; if (d != b * b && !(n & 1)) ans--;
if (ans < 0) ans += ZZQ; cout << ans << endl;
return 0;
}