[TopCoder 12984] TorusSailing（高斯消元主元法优化）

最新推荐文章于 2023-07-24 18:50:40 发布

ixRic

最新推荐文章于 2023-07-24 18:50:40 发布

阅读量580

点赞数

分类专栏： # 高斯消元 # DP 文章标签：高斯消元

本文链接：https://blog.csdn.net/C20190102/article/details/104233900

版权

DP 同时被 2 个专栏收录

31 篇文章 0 订阅

订阅专栏

高斯消元

3 篇文章 0 订阅

订阅专栏

文章目录

题面
分析
代码

题面

Vjudge TorusSailing

分析

首先得到一个DP方程 $dp[i][j]=\dfrac{1}{2}\left(dp[(i + 1)\ \text{mod}\ N][j]+dp[i][(j + 1)\ \text{mod}\ M]\right)+1\quad\quad①$ ~~高消 $\sout{O(N^6)}$ 走人。~~

不能直接高消，我就想按照CF963E Circles of Waiting的做法，由于有很多系数为 $0$ ，那就标记一下不为零的地方，乱搞一下把高消变成 $O(k\cdot n)$ 的。经过超过2小时的尝试，我证明了自己还是太菜太naive。

由于这个环面上可以随意“平移”，不妨把 $(g o a l X, g o a l Y)$ 放到 $(N - 1, M - 1)$ 的地方，那么我们的问题就变成了求 $(N - 1 - g o a l X, M - 1 - g o a l Y)$ 到 $(N - 1, M - 1)$ 的步数期望。

我们把第 $N - 1$ 行和第 $M - 1$ 列的 $d p$ 值作为主元，根据 $①$ 式，我们可以用这些主元由下至上，由右至左表示出其他的 $d p$ 值，然后再由 $①$ 式将第 $0$ 行、第 $0$ 列、第 $N - 1$ 行、第 $M - 1$ 列之间的关系表示成一个 $N + M - 1$ 元 $1$ 次方程组，用高斯消元解这个方程组，时间复杂度 $O(N^3)$ ，进而可以知道所有 $d p$ 值。

代码

点 $(N - 1, M - 1)$ 不要处理，否则会算重。如果访问到 $(N - 1, M - 1)$ 说明起点终点重合，答案是 $0$ 。

#include <algorithm>
#include <cstdio>
#include <cstring>
#include <vector>

const int MAXN = 100;

struct Node {
	double A[MAXN * 2 + 5], B;
}P[MAXN + 5][MAXN + 5];

double Abs(double x) {
	return x > 0 ? x : -x;
}

double A[MAXN * 2 + 5][MAXN * 2 + 5], Res[MAXN * 2 + 5];
	
void Gauss(int N) {
	for (int i = 0; i < N; i++) {
		double tmp = A[i][i];
		for (int j = i; j <= N; j++)
			A[i][j] /= tmp;
		for (int j = i + 1; j < N; j++) {
			tmp = A[j][i];
			for (int k = i; k <= N; k++)
				A[j][k] -= tmp * A[i][k];
		}
	}
	for (int i = N - 1; i >= 0; i--) {
		Res[i] = A[i][N];
		for (int j = i + 1; j < N; j++)
			Res[i] -= A[i][j] * Res[j];
	}
}

class TorusSailing {
	public:
	double expectedTime(int N, int M, int X, int Y) {
		for (int i = 0; i < N - 1; i++)
			P[i][M - 1].A[i] = 1;
		for (int i = 0; i < M - 1; i++)
			P[N - 1][i].A[i + N - 1] = 1;
		for (int i = N - 2; i >= 0; i--)
			for (int j = M - 2; j >= 0; j--) {
				P[i][j].B = 0.5 * (P[i][j + 1].B + P[i + 1][j].B) + 1;
				for (int k = 0; k < N + M - 2; k++)
					P[i][j].A[k] = 0.5 * (P[i][j + 1].A[k] + P[i + 1][j].A[k]);
			}
		for (int i = 0; i < N - 1; i++) {
			A[i][i] = 1, A[i][N + M - 2] = 0.5 * P[i][0].B + 1;
			for (int k = 0; k < N + M - 2; k++)
				A[i][k] -= 0.5 * (P[i][0].A[k] + P[i + 1][M - 1].A[k]);
		}
		for (int i = 0; i < M - 1; i++) {
			int j = i + N - 1;
			A[j][j] = 1, A[j][N + M - 2] = 0.5 * P[0][i].B + 1;
			for (int k = 0; k < N + M - 2; k++)
				A[j][k] -= 0.5 * (P[0][i].A[k] + P[N - 1][i + 1].A[k]);
		}
		Gauss(N + M - 2);
		X = N - 1 - X, Y = M - 1 - Y;
		double Ans = P[X][Y].B;
		for (int i = 0; i < N + M - 2; i++)
			Ans += Res[i] * P[X][Y].A[i];
		return Ans;
	}
};