[JLOI2015] 有意义的字符串 题解

看到这个 $7\times {10}^{18}$ 的模数已经可以摆烂了。

不是，你告诉我这东西跟矩阵快速幂优化 DP 有关系？？

观察到这个题显然不能硬做，因为你显然不能直接算小数部分，而且还有个取模很难受。

所以我们希望把一切的计算都基于整数。

这个时候我们就要思考，有什么东西可以把根号转化为整数。

直接平方再开方显然是不行的，因为你没办法开方取模。

于是，很神奇的东西出现了，你发现这个东西很像求根公式。

即方程 $x^2-b\times x+\frac{b^2-d}{4}=0$ 的两个根分别为 $\frac{b+\sqrt{d}}{2}$, $\frac{b-\sqrt{d}}{2}$。

这个东西显然可以拿来降幂，即 $x^n=b\times x^{n-1}+\frac{d-b^2}{4}\times x^{n-2}$。

所以你就可以直接通过矩阵快速幂求得 $x^n=a\times x^0+b\times x^1$ 的 $a,b$。但是有个问题是 $x^1$ 并不是整数，而前面还有个取了模的系数 $b$ 完全没办法向下取整。

发现一个有趣的事实是，对于这个方程的两个根，他们的和为 $b$。

故我们考虑对于 $f_n=(\frac{b+\sqrt{d}}{2}{)}^n+(\frac{b-\sqrt{d}}{2}{)}^n$ 进行降幂，即 $f_i=b\times f_{i-1}+\frac{d-b^2}{4}\times f_{i-2}$，然后有 $f_0=2,f_1=b$。可以轻松求得 $f_n$。

另外一个有趣的事实是这个奇怪的数据范围 $b^2\le d\le (b+1{)}^2$。可以发现，当 $n$ 为偶数，有 $0\le (\frac{b-\sqrt{d}}{2}{)}^n\le 1$，否则为奇数时，有 $-1\le (\frac{b-\sqrt{d}}{2}{)}^n\le 0$。用人话说就是当 $n$ 为偶数并且 $b^2\ne d$ 时，答案就是 $f_n-1$，否则就是 $f_n$。

#include <bits/stdc++.h>
using namespace std;
#define ll long long
const ll mod = 7528443412579576937ll;
ll add(ll x, ll y){return (1ull * x + 1ull * y) % mod;}
ll mul(ll x, ll y)
{
	ll sum = 0;
	while(y)
	{
		if(y & 1) sum = add(sum, x);
		y >>= 1, x = add(x, x);
	}
	return sum;
}
struct matrix
{
	vector<vector<ll> > a;
	matrix operator *(const matrix &c)const
	{
		matrix now;
		int n = a.size(), m = c.a.size(), l = c.a[0].size();
		now.a.resize(n);
		for (int i = 0; i < n; ++i) now.a[i].resize(l); 
		for (int i = 0; i < n; ++i)
		{
			for (int j = 0; j < l; ++j)
			{
				for (int k = 0; k < m; ++k)
				now.a[i][j] = add(now.a[i][j], mul(a[i][k], c.a[k][j]));
			}
		}
		return now;
	}
	matrix(){}
	matrix(int len1, int len2, int op = 1)
	{
		a.resize(len1);
		for (int i = 0; i < len1; ++i) a[i].resize(len2);
		if(op) for (int i = 0; i < len1; ++i) a[i][i] = 1;
	}
};
ll b, d, n;
matrix ksm(matrix x, ll y)
{
	matrix sum(2, 2);
	while(y)
	{
		if(y & 1) sum = sum * x;
		x = x * x, y >>= 1;
	}
	return sum;
}
int main()
{
	cin >> b >> d >> n;
	if(!n)
	{
		puts("1");
		return 0; 
	}
	matrix now(2, 2), a(1, 2);
	now.a[0][0] = b, now.a[1][0] = (d - b * b) / 4, now.a[0][1] = 1, now.a[1][1] = 0;
	a.a[0][0] = b, a.a[0][1] = 2;
	now = ksm(now, n - 1);
	a = a * now;
	if(b * b != d && n % 2 == 0) a.a[0][0]--;
	cout << a.a[0][0] << endl;
}