Week14 作业 C - Q老师的考验 HDU - 1757

题目描述：

Q老师 对数列有一种非同一般的热爱，尤其是优美的斐波那契数列。

这一天，Q老师 为了增强大家对于斐波那契数列的理解，决定在斐波那契的基础上创建一个新的数列 f(x) 来考一考大家。数列 f(x) 定义如下：

当 x < 10 时，f(x) = x；
当 x ≥ 10 时，f(x) = a0 * f(x-1) + a1 * f(x-2) + a2 * f(x-3) + …… + a9 * f(x-10)，ai 只能为 0 或 1。

Q老师 将给定 a0～a9，以及两个正整数 k m，询问 f(k) % m 的数值大小。

聪明的你能通过 Q老师 的考验吗？（k < 2e9, m < 1e5）

思路：

• 题目给定的是线性递推式，如果使用普通方法求每一个f(x)则复杂度是O(10*K)，无法接受
• 如何优化呢？至少现在我们能想到，f(n)一定能表示成f(0)...f(9)的函数，如果这样不好理解，可以想到我们的高中学过数列，给定数列的初始项，和递推关系，有时候要求其通项公式，如果得到了通项公式，则求f(x)就是O(1)的，那么问题就转换成了如何求f(x)的通项公式
• 其实我们可以把数升级到矩阵，每一次递推就相当于一次矩阵运算，多次递推就相当于迭代矩阵运算，具体的，我们可以推导出如下结果：
• ，矩阵右上角的n-9表示幂次（该图片来自博客https://blog.csdn.net/weixin_44203022/article/details/106330018）

计算矩阵多次幂

• 借助普通快速幂的思想，我们也可以得到矩阵快速幂，只不过一次相乘的复杂度从O(1)变成了O(N^3)，别看O(N^3)很大，主要矛盾在于N不会很大，并且logK也不会很大（logK等于乘法的次数）

• 具体如何实现：封装矩阵类，重载乘法运算符，动态分配矩阵内存（比静态更灵活），快速幂的单位元变成了单位矩阵

代码：

• 之前写的是静态分配数组空间，但是最后期末模拟测试出现了矩阵的大小是由输入的数值决定的，这种情况有两种处理办法：
• 第一种是动态分配空间，如果这样做，不建议使用memset和memcpy，除非你愿意承担意外的风险
• 第二种是开很大的矩阵，但是在计算乘法时，传参N（矩阵的大小），只用一部分空间，参数未必放在重载运算符的函数里，可以用一个类变量来标识，使得重载乘法时使用这个变量
• 下述代码没有内存回收，因为必要性不是很大

#include <cstdio>
#include <iostream>
using namespace std;
long long K, MOD;
//这里定义的矩阵其实是方阵,因为行等于列
struct Matrix
{
	long long **x;
	int N;
	//Constructor
	Matrix(int n)	
	{
		N = n;
		x = new long long* [N + 1];
		for (int i = 1; i <= N; i++)
			x[i] = new long long[N + 1];
		//memset 0
		for (int i = 1; i <= N; i++)
			for (int j = 1; j <= N; j++)
				x[i][j] = 0;
	}

	//可以定义析构函数,但是不是非常有必要

	//operator *
	Matrix operator*(const Matrix& y) const
	{
		Matrix ret(y.N);  //ret==return
		for (int i = 1; i <= N; ++i)
			for (int j = 1; j <= N; ++j)
			{
				ret.x[i][j] = 0;
				for (int k = 1; k <= N; ++k)
					ret.x[i][j] = (ret.x[i][j] + (x[i][k] % MOD) * (y.x[k][j] % MOD) % MOD) % MOD;

			}
		return ret;
	}

	//copy Constructor
	//不建议使用memcpy,因为是动态申请的,对指针memset或memcpy可能出问题
	Matrix(const Matrix& t)  
	{
		N = t.N;
		x = new long long* [N + 1];
		for (int i = 1; i <= N; i++)
			x[i] = new long long[N + 1];
		for (int i = 1; i <= N; i++)
			for (int j = 1; j <= N; j++)
				x[i][j] = t.x[i][j];
	}
};
Matrix mqpow(Matrix a, int x)
{
	Matrix ret(a.N);
	//置为单位矩阵
	for (int i = 1; i <= ret.N; i++)
		ret.x[i][i] = 1;
	//迭代矩阵快速幂
	while (x)
	{
		if (x & 1) ret = ret * a;
		a = a * a;
		x >>= 1;
	}
	return ret;
}
int main()
{
	while (cin >> K >> MOD)
	{
		Matrix A(10);
		for (int i = 1; i <= 10; i++)
			cin >> A.x[1][i];
		for (int i = 2; i <= 10; i++)
			A.x[i][i - 1] = 1;
		Matrix HalfAns = mqpow(A, K - 9);
		long long FinalAns = 0;
		for (int i = 1; i <= 10; i++)
			FinalAns = (FinalAns + HalfAns.x[1][i] * (10LL - i) % MOD) % MOD;
		cout << FinalAns << endl;
	}
	return 0;
}