自动搬运

查看原文

来自洛谷，原作者为

	Helenty 江风引雨 / 恣睢浇灭少年抱有的幻想 / 断弦声里 / 无计唤停没有结果的爆搜 / 此去经年 / 是否还能放下封存的回忆

搬运于2025-08-24 22:25:02，当前版本为作者最后更新于2025-04-23 10:28:17，作者可能在搬运后再次修改，您可在原文处查看最新版

自动搬运只会搬运当前题目点赞数最高的题解，您可前往洛谷题解查看更多

以下是正文

显然的高斯消元，但要优化。

消元时可以只消第 $i$ 行上下 $m$ 行，消元时也只消有元的那 $m$ 列，复杂度是 $O(nm^3)$。

注意消元时如果使用高斯约旦消元的话，那么会破坏矩阵性质，要用有会带的那种消元。

在记录的时候可以把矩阵 $i$ 行 $j$ 列记到 $i,j-i+m$ 这样子就只用记录 $2nm$ 个值。

在矩阵消元时如果遇到 $0$，当且仅当这个点和第一行不连通，可以 continue 。回代时也是，输出时可能要特判。

代码如下，轻松最优解前三。

#include <bits/stdc++.h>
using namespace std;

const int N = 21, M = 10005;
const double eps = 1e-9; 

int n, m, c;              
char s[M][N];             
double u, d, l, r;     
double mp[N*M][N*4];      
double ans[N*M];    
double p; 

int main() {

    scanf("%d%d%lf%lf%lf%lf", &m, &n, &u, &d, &l, &r);
    
    for (int i = 1; i <= n; i++)
        scanf("%s", s[i] + 1);
    
    for (int i = 1; i <= n; i++) {
        for (int j = 1; j <= m; j++) {
            if (s[i][j] == 'X') { 
                mp[(i-1)*m + j][m] = 1;
            } else {
                if (i == 1) c++; 
                
                double g = 0; 
                if (s[i][j] == '.') {
                    if (i != n && s[i+1][j] != 'X') g += d;
                    if (i != 1 && s[i-1][j] != 'X') g += u;
                    if (j != 1 && s[i][j-1] != 'X') g += l;
                    if (j != m && s[i][j+1] != 'X') g += r;
                    
                    if (i != n && s[i+1][j] != 'X')
                        mp[(i)*m + j][0] = d / g;
                    if (i != 1 && s[i-1][j] != 'X')
                        mp[(i-2)*m + j][2*m] = u / g;
                    if (j != 1 && s[i][j-1] != 'X')
                        mp[(i-1)*m + j-1][m+1] = l / g;
                    if (j != m && s[i][j+1] != 'X')
                        mp[(i-1)*m + j+1][m-1] = r / g;
                }
                
                mp[(i-1)*m + j][m] = -1;
            }
        }
    }
    
    for (int i = 1; i <= m; i++)
        if (s[1][i] == '.')
            ans[i] = -1.0 / c;
    
    // 带状高斯消元（前向消元）
    for (int i = 1; i <= n*m; i++) {
        if (fabs(mp[i][m]) <= eps)
            continue; // 跳过零行
        
        // 消去下方m行内的元素
        for (int j = i + 1; j <= min(i + m, n*m); j++) {
            double ratio = mp[j][m - j + i] / mp[i][m];
            
            // 更新系数矩阵
            for (int k = m; k <= 2*m; k++)
                mp[j][k - j + i] -= ratio * mp[i][k];
            
            // 更新解向量
            ans[j] -= ans[i] * ratio;
        }
    }
    
    // 回代求解
    for (int i = n*m; i >= 1; i--) {
        for (int j = i + 1; j <= min(i + m, n*m); j++)
            ans[i] -= mp[i][m + j - i] * ans[j];
        
        if (fabs(mp[i][m]) > eps)
            ans[i] /= mp[i][m]; // 归一化
        else
            ans[i] = 0;         // 零行处理
    }
    
    // 输出结果并校验概率总和
    p = 0;
    for (int i = 1; i <= n; i++)
        for (int j = 1; j <= m; j++)
            if (s[i][j] == 'T')
                printf("%.10lf\n", ans[(i-1)*m + j]), p += ans[(i-1)*m + j];
    
    // 概率总和校验
    assert(fabs(1 - p) <= eps);
    
    return 0;
}